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The Influence of Biaxial Mean Stresses on 
the Failure of Tubes by Fatigue 


H. L. COX, M.A., F.R.Ae.S., A-M.I-Mech.E. and N. B. OWEN, A.M.1.Mech.E. 


(National Physical Laboratory) 


SUMMARY: Thin-walled tubes, § in. in diameter, of three hard aluminium alloys and 
of mild steel have been tested in fatigue under three systems of alternating stresses 
while subjected to biaxial mean tensions imposed by means of internal pressure. 
In fatigue under direct or bending stresses the hoop tension in the walls of the 
aluminium alloy tubes did not seriously reduce the fatigue endurance, but it did 
markedly affect the mode and rate of crack propagation; cracks initially transverse 
to the tube axis tended to develop very rapidly in the axial direction. This tendency 
was present under both fluid and gas pressure, and under gas pressure the cracks 
propagated so fast that the test piece was often blown to pieces before the gas 
pressure fell by leakage through the cracks. The gradual taper in wall thickness 
along the fillets joining the test section to the enlarged ends offered no barrier to 
propagation of the axial cracks and the whole test piece, including its enlarged ends, 
was often shattered. Propagation of the axial cracks was preventible by sufficiently 
reducing the fillet radius, or by a ring glued on. Under alternating torsion, both 
endurance and mode of failure were affected by internal gas pressure. The initial 
fatigue crack, either circumferential or axial, often extended over a length com- 
parable with the diameter of the tube, except under low ranges of shear stress when 
the crack length was sometimes very short. At each end the cracks forked in a 
characteristic manner and under moderate gas pressure the portions of the wall 
between the prongs of the fork were blown outwards. Under high pressure 
explosive failure and fragmentation often occurred. Mild steel under alternating 
torsion with internal gas pressure exhibited the same modes of failure, and two or 
more fatigue cracks were often formed simultaneously. Under high hoop tension, 
cracks propagated rapidly and one test piece, after two million cycles endurance, 
failed by exploding. However, no mild steel test piece was fragmented. A tentative 
explanation is offered of the reason why, in tubes of small diameter, rather short 
fatigue cracks may be expected to lead to fast fracture under the static loading. 
Attention is drawn to the inference that the initial fatigue crack itself must develop 
very quickly to considerable length. 


1. Introduction 


An investigation has been made into the effect of biaxial mean stress on failure 
by fatigue. 


Test pieces were made in the form of tubes, of § in. diameter and 0-015 in. wall 
thickness, and the first tests were made under static hoop tension, imposed by internal 
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pressure, and pulsating axial tension. The axial loading was imposed by means of the 
N.P.L. 1500 Ib. Slipping Clutch Fatigue Machine, the mean loading being set to half 
the range imposed, diminished by the axial loading due to the internal pressure. For 
convenience the internal pressure was obtained from a cylinder of gas. 


In the first test under these conditions the test piece of aluminium alloy to 
D.T.D. 683 endured 27,000 cycles and then exploded, fragments being blown all over 
the place (VCF 3, Fig. 10 on page 23). This was surprising, because the precaution 
of closing the valve of the gas cylinder had been taken and the energy available to 
fragment the test section was only that of the gas within the tube under test. 


Two tines developed out of this first test; one was to investigate the conditions 
under which explosive failure occurred, and the other to ascertain the degree in which 
this type of failure was ascribable to the use of gas as the pressure medium. 


In prosecution of the second line, tests were made using oil as the pressure 
medium and other tests with a reduced volume of gas. The three conclusions from 
these tests were quite clear; first, that fragmentation was due to the energy in the gas 
and that, by sufficient reduction of the gas volume, failures similar to those resulting 
when oil was used as the pressure medium could be obtained; second, that the 
tendency to fragment was inherent in the material itself, being far greater in D.T.D. 
683 than in B.S. L 65; and third, that the tendency to fragment arose from a liability 
of the materials to develop several axial cracks, which often extended well beyond the 
thin-walled section and even right through the much thicker enlarged ends. 


The gradual thickening of the wall associated with a generous radius of transition 
fillet appears to offer little resistance to the propagation of these axial cracks, which 
can penetrate into sections where the wall thickness is ten times greater, whereas a 
more abrupt change of section, afforded by a small fillet radius, does prevent the crack 
spreading even into a section only four times thicker. A very small fillet radius, or a 
separate ring glued on, will result in circumferential cracking, after a shorter endurance 
but without catastrophic failure. 


In respect of the comparison between failures under gas and oil pressure media, 
our conclusions have since been fully substantiated elsewhere in tests on much larger 
cylinders. It is notable that the tendency to fragment depends on the total energy in 
the gas and this energy is related not to the gas pressure p itself but to the hoop 
tension, pD/(2t), to which that pressure gives rise. In short, the complete fragmen- 
tation of our cylinders with D/t of about 40 under 1000 Ib./in.2 pressure may presage 
similar fragmentation of cylinders with D/t of about 3000 under a pressure of 13 
Ib./in.2 


The second and third conclusions led us, in the first and main line of investigation, 
to compare the behaviour of the three aluminium alloys, D.T.D. 683, B.S. L 65, B.S. 
L 64, and a mild steel under several systems of loading. The general result of the tests 
has been to confirm marked differences between the materials in their liability to 
fragment; D.T.D. 683 is by far the worst, B.S. L 65 is markedly less liable, but 
apparently slightly worse than B.S. L 64, while mild steel has not yet been made to 
fragment. 
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At the same time it soon became apparent that explosive failure as such was not 
inevitably associated with the tendency to fragment. Explosive failure seemed to 
depend not entirely on metallurgical factors but also very markedly on the mode of 
stressing. In short, we were forced to conclude that, under suitable conditions, even 
mild steel might fail explosively and we have indeed succeeded in causing mild steel 
test pieces to explode after endurances ranging from 10‘ up to 2 million cycles. 


It appears that the type of stressing most likely to lead to explosive failure is that 
in which planes subjected to maximum range of shear stress are also subject to maxi- 
mum static direct stress; this observation might have been anticipated—though it was 
not—from previous experience. The influence of metallurgical factors inherent in 
extruded or rolled bars is also marked and this influence may sometimes predominate. 


An empirical explanation of explosive failure of tubes is offered, which appears 
to correlate the observations both in the present tests and in many tests on larger 
cylinders made elsewhere. However, this explanation merely relates the occurrence 
of fast failure to the dimensions of the tube and indicates how the unstable crack 
length may be related to the tube dimensions. It leaves open two other questions; 
one, what determines the unstable crack length in a flat plate? and the other, what 
determines the growth of the fatigue crack up to the unstable length? 


Both questions are fundamental but, while the first is the subject of much intensive 
study at present in progress, the second is not receiving so much attention. The results 
of the tests described in this paper clearly indicate that the propagation of the fatigue 
crack up to the unstable length is often extremely fast. Indeed, it would scarcely be 
too much to claim that the initial fatigue crack may sometimes spring into existence 
at a length approaching or exceeding the unstable crack length. It is this occurrence 
which results in explosive failure and fragmentation. 


2. Materials Tested 


Tests were made on five aluminium alloys, supplied in | in. diameter extruded 
bars, two to specification D.T.D. 683, two to B.S. L 65 and one to B.S. L 64. The 
chemical compositions of these alloys and their heat treatments are given in 
Table I(a). 

Tests were also made on a mild steel of 30 tons/in.* tensile strength, also in the 
form of | in. diameter bar, the chemical composition of which is given in Table I (6). 

The mechanical properties of all the materials used are given in Table II: the 
0-1 per cent proof stress is quoted for the aluminium alloys and the yield stress for the 


steel. Table II also includes the tensile hoop stresses obtained from static bursting 
tests described in Section 4. 


3. Preparation of Test Pieces 

The tubular test pieces used for the main series of tests on aluminium alloys are 
shown in Figs. 1 (a) and 1 (6). They were made by first drilling an undersize hole, 
which was then finished to size by boring; in order to achieve the desired high quality 
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FIGURE |. The test pieces used (all dimensions are in inches). 


of finish in these deep holes, extremely fine cuts were taken in the boring process. 
For the piece shown in Fig. | (a) the threaded ends were then formed and finally the 
test portion and fillets were turned coaxial with the bore; for the piece shown in Fig. 
1 (6) the test portion and fillets were turned coaxial with the bore, before the square 
ends were machined. The steel pieces of this size were similar, but the wall thickness 
was 0-010 in. 


To have polished the bore would have meant too great a risk of leaving it out of 
round; therefore, to maintain uniform wall thickness, both inner and outer surfaces 
were left as machined. The machine finish was, in fact, fairly good, and the influence 
of the surface irregularities on the fatigue strength was unlikely to have been great. 
In any case, this influence should have favoured the development of circumferential 
cracks, because all the turning marks ran in that direction, whereas in many tests axial 
cracks predominated. 


Figure 1 (d) shows the form of a smaller test piece used in a later series of tests on 
the steel. In this case, since the required wall thickness was only 0-005 in., the test 
piece was redesigned as an open tube and the bore was lapped to allow precise location 
on the mandrel used to support it during the machining of the outer form. The test 
portion was polished circumferentially with 00 emery paper. 
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4. Bursting by Internal Pressure Only 


To afford a comparison with the fatigue tests described later on, a series of tests 
was done in which test pieces were burst by the application of internal gas pressure 
only. The test piece was first filled with gas from a bottle at 2000 Ib./in.? pressure; 
then the gas supply was disconnected and oil pressure was used to compress the gas 
further until failure occurred. The gas volume involved in these tests was roughly the 
same as in subsequent tests under fatigue loading. 


Tests were made on § in. diameter tubes in the materials D.T.D. 683, B.S. L 65 
and B.S. L 64; a4 in. tube of steel was also tested under gas pressure, and a § in. tube 
of D.T.D. 683 was burst by oil pressure. 


Photographs of the broken test pieces are shown in Fig. 8 on page 21. The values 
of ultimate tensile hoop stresses obtained are included in Table II; the value quoted 
for material VCF is the mean of three tests made with both the gas and oil pressure 
media; it is notable that the difference between the greatest and least of the three 
results was only 2 tons/in.2._ The values quoted for the other materials are based on a 
single test made with a gas pressure medium. 


5. Fatigue Tests Under Biaxial Mean Stresses 
5.1 MODE OF TEST 


The mode of test adopted was to subject a tubular test piece to static internal 
pressure and to superimpose a fatigue loading of either direct stress, bending or 
torsion. The N.P.L. 1500 Ib. Slipping Clutch Fatigue Machine was used for the direct 
stress fatigue tests and the N.P.L. Combined Stress Fatigue Machines, both large and 
small, were used for the torsion and bending fatigue tests. 


The dimensions of the test pieces shown in Fig. 1 were determined by the 
capacities of the fatigue testing machines available and by the highest pressure (2000 
lb./in.2) readily available from gas cylinders. 


To reduce the number of joints to be made gas tight, the test piece was made as a 
tube with a closed end. For the direct stress tests, access for the pressure medium to the 
bore of the tube was provided through the lower machine adaptor; in the torsion and 
bending tests, the supply pipe was connected straight to the test piece, through a 
suitable screwed connector, and this end of the test piece was mounted at the frame 
end of the machine. In the early tests a soft copper washer was used as a seal at the 
test piece end, but this was later replaced by a seal of hemp and jointing compound, 
which proved more satisfactory. The gas cylinder was connected to the test piece 
through a copper pipe, of 0-125 in. bore and about 7 ft. long, and the pressure applied 
was read by a gauge at the pipe end immediately adjacent to the gas cylinder. After 
the desired pressure had been attained, the valve on the cylinder was almost closed 
again, leaving just a sufficient opening to annul the very slight leakage in the system. 
Consequently, when the test piece failed the pressure fell at once to zero, and the gas 
released was restricted to that in the test piece and connecting pipe. 
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The procedure adopted for testing was to mount a test piece in the machine, 
test the system for pressure tightness, and then to set the pressure to the required value. 
When the fatigue loading was pulsating direct stress, the mean axial load was applied, 
this load being adjusted to the value which, together with the axial load due to the 
internal pressure, would afford the desired mean axial stress: the machine was then 
set running at the desired load range. In the tests in the combined stress machine, the 
gas pressure set up both hoop tension and axial tension half as great. When the 
pressure had been set to the desired value and maintained for a few minutes, the 
machine was set running at the desired range of load. 


Examination of the fractured surfaces of some of the early pieces tested in direct 
stress indicated that these surfaces had been forced together by the sudden reduction 
of mean load resulting from the release of gas pressure when the test piece broke. 
This effect was overcome by fitting stops to limit the axial motion of the machine 
loading frame. 


5.2. Direct STRESS PULSATING TESTS 


The tests made under pulsating direct stress fatigue loading, with a superimposed 
hoop tensile stress, were sub-divided into two main groups. One series, to determine 
the effect of the hoop stress on the mode of failure, was run at a constant value of 
pulsating stress with the pieces loaded at several values of static hoop stress; the other 
series, to determine the effect of the pulsating stress, was run at a low value of hoop 
stress over a range of pulsating stresses. 


The results of these tests are given in Tables III and IV, and the S-N curves are 
plotted in Figs. 2, 3, 4 and 5. Only a few check tests were made on the steel and these 
results are not shown plotted but are given in Table IV only. Photographs of the 
broken test pieces in all the materials tested are given in Fig. 10 on page 23. 


5.3. TESTS WITH VARYING ENERGY CONTENT OF PRESSURE MEDIUM 


The pulsating direct stress tests mentioned in Section 5.2 were made using a gas 
pressure medium; the energy released at fracture by expansion of the gas within the 
test piece was of the order of 100 ft.lb. Clearly this release of energy would be 
expected to enhance the tendency to shatter and, accordingly, to determine how far 
shatter depended upon this release of energy, a number of tests were made with a 
pressure medium of varying energy content. 

Tests were made on three pieces of D.T.D. 683, to Fig.1 (a), in which the internal 
volume was reduced by fitting plugs of varying clearance to the bore of the test pieces. 
These reduced the volume to about 20 per cent, 3 per cent and 0-3 per cent of the total, 
while leaving the hoop stress unaltered; the results of these tests are given in Table V. 
Three further test pieces of D.T.D. 683 were tested full of oil and, also, one piece of 
B.S. L 65; the results of these tests are included in Table III, and are plotted in Fig. 2 
for comparison. Photographs of these broken pieces, which show the effect of energy 
content on the mode of failure, are given in Fig. 9 on page 22. 
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Ficure 5. Material: B.S. L 64. N.P.L. test mark: VCH. 
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5.4. EFFECT OF STRESS GRADIENT 


A limited number of tests were made on material D.T.D. 683, to check the effect 
of the transition radii at the ends of the test portion on the propagation of the 
longitudinal cracks. From the other tests it was known that these cracks developed 
extremely fast in D.T.D. 683, particularly at high values of pulsating direct stress and 
hoop stress, and propagated well into the enlarged ends of the test pieces, as shown in 
Fig. 10 on page 23. The form of test piece used in this series of tests is shown in 
Fig. 1(c) and tests were made with fillet radii of 4, 3, } and 14 in.; pieces with 4 in. 
radii had already been covered in the other tests. In addition, a test was made on a 
piece to Fig. 1(a), but with a ring { in. wide Redux-bonded around the centre of the 
test portion. 


The pieces were all tested at about the same condition of loading, a high range of 
pulsating direct stress and a high hoop stress. The results of the tests are given in 
Table VII and the endurances to fracture are plotted in Fig. 6; included in this figure, 
but not in the table, are the results of some earlier tests made with varying fillet radii. 
Photographs of the broken test pieces are shown in Fig. 14 on page 27, and Fig. 7 
shows a diagram of the nominal stress gradients associated with the fillet radii tested. 


5.5. ALTERNATING BENDING TESTS 


A few check tests were made under alternating bending with superimposed tensile 
hoop stress and axial tension on pieces to Fig. 1 (6) in all the aluminium alloys. The 
results, given in Table VI, are plotted as S-N curves in Figs. 3,4 and 5. Photographs 
of the broken pieces are shown in Fig. 11 on page 24: the mode of failure of the 
piece VCH 8 (B.S. L 64) relates to a test made under a low hoop stress, but another 
test, made under a hoop stress of 45 per cent, resulted in a fracture almost 
indistinguishable from that shown. 
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FiGure 6. Effect of stress gradient on crack penetration. D.T.D. 683. 
@ N.P.L. test mark VCF 
+ N.P.L. test mark WMS 
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Ficure 7. Effect of stress gradient on crack propagation under fatigue. 
(x is the distance in inches from the parallel section.) 


For comparison the results of rotating-bending fatigue tests on materials D.T.D. 
683, B.S. L 64 and B.S. L 65 are given in Table VIII and are also included in Figs. 
3,4 and 5: the results for the mild steel are given in Table VIII only. 


5.6. ALTERNATING TORSION TESTS 


A number of tests under alternating torsional stress with superimposed tensile 
hoop stress and axial tension were made on all materials; the pieces were to Fig. 1 (0). 
The first tests were made with a constant alternating stress at different values of static 
stresses. Because the results obtained were most unexpected, additional tests were 
made on D.T.D. 683 and B.S. L 65 at various levels of alternating torsion, with and 
without static stresses. The results of all these tests are given in Table IX, and the S-N 
curves are plotted in Figs. 3, 4 and 5; Figs. 12, 13 and 16 show photographs 
of some of these 3 in. diameter broken test pieces. 


6. Unstable Crack Length 


In Appendix I an empirical theory is developed which indicates how the unstable 
crack length may be expected to diminish as the size of the tube is reduced. The values 
of unstable crack length indicated by this theory are in fair agreement with observ- 
ations made on larger tubes elsewhere (see Appendix I). 

Observation of the crack lengths developed in % in. tubes of mild steel under 
fatigue loading (piece BNZ 13F5, Fig. 16 on page 28) suggested that in smaller tubes 
the initial fatigue crack might prove unstable. In order to test this inference, a 
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number of mild steel test pieces of the form shown in Fig. I(d), 0-25 in. bore and 
(0-005 in. wall thickness, were tested under alternating torsion and internal pressure 
affording hoop and axial static tensions. 


The results of the tests are given in Table X(a) and a photograph of the pieces 
after test is given in Fig. 15. Included in Table X(a) are the apparent lengths of the 
axial fatigue cracks as measured from the broken pieces, using a travelling micro- 
scope, and the unstable crack lengths as computed according to Appendix I, using 
the experimental value of the hoop tensile stress, namely 28:5 tons /in.* 


Because of the uncertainty of the stress combination that these particular pieces 
would withstand, the stress applied to the first piece was increased in stages, the hoop 
stress being first raised to the maximum gas pressure available before the torsional 
stress was increased. The stresses quoted are nominal values based on the mean wall 
thickness of the tube. 


Subsequently a series of pieces of mild steel of the larger size (Fig. 1 (5), but 
having the wall thickness 0-010 in.) were tested under alternating torsion, with super- 
imposed static stresses set up by internal gas pressure. The results of these tests are 
given in Table X(b) and a photograph of some of the pieces after test is given in 
Fig. 16. 

7. Conclusion 


The evidence afforded by the present tests has to be reviewed in relation to two 
separate questions :— 


(a) In what way and to what extent does biaxial mean stress affect failure by 
fatigue ? 


(b) Under what conditions does fatigue failure develop almost instantaneously 
into complete rupture? 


For convenience of comparison in relation to these questions, the hoop tensile 
and fatigue stresses quoted with the photographs (Figs. 10 to 16) are there expressed 
as ratios to the tensile strengths of the materials. The actual values are given in the 
tables. 


Under direct stress or bending, high hoop stress due to internal pressure does not 
materially reduce endurance under high fatigue stresses; nor is the endurance affected 
by the energy content of the pressure medium. Under torsion, however, internal 
pressure does markedly reduce endurance; for instance, one piece (VCF 29) after 
an endurance of 10 million cycles without internal pressure, failed when internal 
pressure was applied after a further 350,000 cycles of the same range of shear stress. 


Sufficiently high internal pressure may markedly affect the mode of failure by 
fatigue and cause it to develop almost instantaneously into complete rupture. 


Under static loading, fast fracture is known to develop once a crack or slit 
exceeds a certain unstable length. Under any particular loading condition the unstable 
length is a characteristic of the material; but, for tubes of any one material, the length 
is shorter the smaller the diameter of the tube and the higher the hoop tensile stress. 
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D.T.D. 683 (VCF 23) 
Burst by gas pressure 
D.T.D. 683 (VCF 27) 
Burst by oil pressure 


B.S. L 65 (MAK 15) B.S. L 64 (VCH 13) 
Burst by gas pressure Burst by gas pressure 


FiGurRe 8. Cracks due to static pressure. 
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D.T.D. 683 (VCF 19) D.T.D. 683 (VCF 20) 


Gas pressure 
18°, of total volume 3 


Gas pressure 
of total volume 


D.T.D. 683 (VCF 6) 
Oil pressure 


D.T.D. 683 (VCF 22) 
Gas pressure 
of total volume 


B.S. L. 65 (VCG 7) 
Oil pressure 


Ficure 9. Effect of energy content of pressure medium on crack propagation under fatigue. 
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D.T.D. 683 (VCF 
Hoop stress 10 
Range of direct stress 65 


D.T.D. 683 (VCF 3) 
Hoop stress 25 
Range of direct stress 65 


Right: 
B.S. L 64 (VCH 1) 
Hoop stress 45 
Range of direct stress 70 


B.S. L 65 (VCG 6) Mild steel (BNZ 13 F 1) 
Hoop stress 55 Hoop stress 35 
Range of direct stress 85 Range of direct stress 60 


FiGuret 10. Fatigue cracks under pulsating direct stress and static hoop tension. 
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D.T.D. 683 (VCF 14) 
Hoop stress 40 
Range of bending stress 75 


D.T.D. 683 (VCF 15) B.S.L 65 (VCG 15) B.S. L 64 (VCH 8) 
Hoop stress 25 Hoop stress 50 Hoop stress 15 
Range of bending stress 80 Range of bending stress 95 Range of bending stress 85 


Ficure Il. Fatigue cracks under alternating bending and static hoop tension. 
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Left 
D.T.D. 683 (VCF 18) 


Hoop stress 10 
Range of torsion stress 40 
Right 
B.S.L. 65 (VCG 16) 


Hoop stress 35 


Range of torsion stress 50 


D.T.D. 683 (VCF 13) B.S. L 65 (VCG 13) 
Hoop stress 40 Hoop stress 50 
Range of torsion stress 45% 


Range of torsion stress 34 


Left 
D.T.D. 683 
(VCF 16) 


1, Hoop stress 


25 


Range of torsion stress 40 


Right 
D.T.D. 683 
(VCF 34) 


Hoop stress 4§ 
Range of torsion stress 10 


FiGURE 12. Fatigue cracks under alternating torsion and static hoop tension 
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B.S. L 65 B.S. L 64 B.S. L 64 
(VCG 17) (VCH 7) (VCH 12) 
Hoop stress 10 Hoop stress 15 Hoop stress 15 


Range of torsion stress 45 Range of torsion stress 45 Range of torsion stress 45 


BS: L. 65 B.S. L 64 B.S. L 64 
(VCG 14) (VCH 11) (VCH 10) 
Hoop stress 15 Hoop stress 30 Hoop stress 50 
Range of torsion Range of torsion stress 45 Range of torsion stress 45 
stress 45 
FIGURE 13. 


Fatigue cracks under alternating torsion and static hoop tension. 
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Fictre 1S, Fatigue cracks under alternating torsion and static hoop tension in mild steel 
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B.S. L 65 B.S. L 64 B.S. L 64 

(VCG 17) (VCH 7) (VCH 12) 

Hoop stress 10 Hoop stress 15 Hoop stress 15 
Range of torsion stress 45 Range of torsion stress 45 Range of torsion stress 45 


B.S. L 65 

(VCG 14) 
Hoop stress 15 
Range of torsion 


stress 


FIGURE 13. 


B.S. L 64 B.S. L 64 

(VCH 11) (VCH 10) 

Hoop stress 30 Hoop stress 50 
Range of torsion stress 45 Range of torsion stress 45 


Fatigue cracks under alternating torsion and static hoop tension. 
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Fatigue cracks under alternating torsion and static hoop tension in mild steel. 


Range of torsion 


stress 85 
FIGURE 15 
1961 27 


February 


FATIGU FAILURI OF} TUBES 
‘ iz 
> 


4: AND NN. OWEN 


Right: 
BNE 2F 22 
No hoop stress 


Range of torsion stress 75 


Lefi 

BNZ2F 17 

No hoop stress 

Range of torsion stress 65 


BNZ 2 F 21 
Hoop stress 65 
Range of torsion stress 45 

ue BNZ 2 F 13 
Hoop stress 65 
Range of torsion stress 55 


Right: 


L is the 
total length 


of the 
axial crack 


BNZ 2 F 23 BNZ 13 F 5 


Hoop stress 70 


Hoop stress 50 
Range of torsion stress 55 sahara 
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Ficure 16. Fatigue cracks under alternating torsion in mild steel with and without static 


hoop tension. 
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Accordingly, in tubes under internal pressure, the first fatigue crack, if it develops 
rapidly enough, may lead at once to fast fracture under the static loading. In con- 
sequence, all the complete failures in the present tests, and in particular the explosive 
ones, look not unlike failures under internal pressure alone. 


However, complete rupture can result under fatigue loading only if the initial 
fatigue crack develops to a length greater than the unstable length appropriate to the 
material, the size of the tube and the mean loading applied, before the gas pressure 
falls to zero. The results of the present tests demonstrate that sufficiently rapid 
development is possible in both hard and soft materials; they demonstrate further, 
on the one hand, that the unstable crack length may be quite short, and, on the other, 
that development of the fatigue crack to this length must be extremely rapid. 


The present tests are too few to justify strictly quantitative assessment of these 
two factors but they do afford some further evidence of a qualitative nature with 
regard to the effects of surface finish and other metallurgical factors. 


Under the conditions of the torsion tests (the hoop tensile stress being twice the 
axial mean tensile stress), fast fracture would be expected to develop more readily 
from axial than from circumferential cracks. The finish of the tubular test pieces 
certainly favoured the latter mode and, in the hard aluminium alloys, the first crack 
was often circumferential. However, when the internal pressure was high, these 
circumferential cracks often branched and developed into axial ones; comparison 
between D.T.D. 683, B.S. L 65 and B.S. L 64 in this respect strongly suggests that 
metallurgical factors, perhaps associated with the extrusion process, may be respon- 
sible for this tendency. In the mild steel pieces of § in. diameter, circumferential 
cracks developed in all tests made without a hoop stress, whereas with a hoop stress 
the cracks were mainly axial. 


Once an axial crack has got going, only an abrupt change of section will stop 
it spreading, but a glued-on ring is equally as effective as a small radius fillet. The 
concentration of stress at such abrupt change of section may reduce the fatigue 
endurance, but the resulting crack will be circumferential and, accordingly, far less 
likely to lead to catastrophic failure. 


It remains to comment on the results of the tests on the material to D.T.D. 683 
under alternating torsion, which are particularly disquieting. Even without mean 
stresses due to internal pressure, the fatigue resistance to this form of loading was 
much lower than the estimate based on the fatigue resistance under direct stress or in 
bending; when biaxial mean stresses were imposed, fatigue failures resulted under 
ranges of shear stress as low as +1°8 tons/in.’ 

The fatigue resistance of material to B.S. L 65 under alternating torsion is also 
low by comparison with its resistance to direct or bending stresses but, in this material, 
biaxial mean stress did not materially affect the strength at long endurances, although 
the endurances at high stresses were appreciably reduced. 

From the literature it certainly appears that the fatigue resistance of the hard 
aluminium alloys under alternating torsion is often low by comparison with their 
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resistance under alternating bending, ratios as low as 0-4 having been reported; but 
the reasons for this weakness seem not to have been explored. 


To supplement the evidence presented here, further tests under alternating 
torsion have been made on two further samples of material to D.T.D. 683 in tubular 
form and on solid pieces of all three materials. The results so far indicate that, in 
tubular form, all three sources of material are equally poor and that their low strength 
is still further reduced by biaxial stresses imposed by internal pressure. One of the 
materials also gave much the same results when tested in the form of solid bar, 
whereas the material VCF and the third material in solid bar form exhibited under 
alternating torsion a resistance fully concordant with their resistance under alternating 
bending. 


The differences in behaviour do not appear to be due to differences in surface 
finish, but it is possible that they may be due to internal stresses. However, this 
explanation seems likely to prove satisfactory only if the high fatigue resistance in 
solid bar form is attributable to compressive mean stresses in the surface layers. 
Accordingly, it appears probable on the evidence at present available that the fatigue 
resistance of the hard aluminium alloys under alternating shear stresses may prove 
abnormally low when the planes subjected to maximum shear are also subjected to 
Static tensile stresses. 
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Appendix I 


UNSTABLE CRACK LENGTH IN A TUBE 


1. Increase of Hoop Tension at the Ends of an Axial Slit 


When an axial slit of length a is made in the wall of a tube, of radius R and wall 
thickness f, which is under internal pressure p, the circumferential tension, originally 
f,=p(R/t) everywhere, is reduced to zero at the edges of the slit. It follows that the 
mean circumferential tension in the way of the slit is reduced by Af,, where A is some 
ratio between 0 and | which must tend to zero with the ratio a/R. 


By analysis of the diffusion of loading applied along the edges of the slit outlined in 
Appendix II, it appears reasonable to assume A=u(a/R), where the factor u is 
about 1/8. 


The internal pressure being everywhere maintained at the value p, over the segment 
of the tube in the way of the slit, a part » (a/R) p of the total pressure is then unbalanced 
by circumferential tension and each of the sections of the tube through the ends of the slit 
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are subjected to shear loading F=4n(a/R)pa and bending moment M= 
1/{12u (a/R) per unit length distributed uniformly* round the perimeter. These 
loadings give rise to an extra circumferential tension f, in these sections where (see 
Appendix III) 


f,=[2-6F (Rt)-'/2+3-3M (Rt)-') (R/0). 
Substituting for F and M, and writing f, for pR/t, 


(1-3+0-36), : 
where 6=a/(Rt)'’? 


2. Effect of Stiffening and of the Periodic Components of Stresses Due to the Slit 

If the tube is stiffened by sturdy rings closely adjacent to the ends of the slit, both 
shear loading and bending moment may be virtually annulled and f, may be reduced 
almost to zero. Similarly, closely spaced sturdy stringers would tend to reduce the 
bending moment but, in practical constructions, the influence of the stringers is likely to 
be far less than the effect of the rings. On the other hand, the effect of the stringers may 
be more or less independent of the position of the slit, whereas the rings will have but 
slight effect unless the ends of the slit are fairly close to them, say, within a distance 
2 (Rr)! 


The mean radial shear loading F and mean bending moment M are, of course, only 
the first (constant) terms in Fourier series which would represent the total loading on the 
outer parts of the tube arising from the disturbance of stress in the way of the slit; and, 
of course, the periodic terms in these series also give rise to stresses in the outer parts. 
However, the heart of plate stresses arising from the periodic terms are mainly shear 
stresses on the axial and circumferential sections and these terms do not contribute at all 
to the mean hoop tension (or rather their contribution is only second order). Moreover, 
the bending strains in the tube walls arising from the periodic terms are limited by their 
association with shear strains in the plane of the wall, and the bending stresses will be 
further limited by the stiffeners. In short, the stresses arising from the periodic terms 
are not likely to be very great and, in any case, these periodic terms do not affect 
the value of f,,. 


3. Interpretation and Application of Formula (\) 


Formula (1) represents the increase of nominal circumferential tension at the ends of 
the slit due simply to the curvature of the tube; if R is infinite, f, is zero. Whereas the 
nominal tension at the ends of the slit in a flat plate would be f,, in a tube it is f,+f,. 


The actual local concentration of stress at the ends of the slit is probably much the 
same in both flat plate and tube, at least if R/t is of order 100. Therefore, in comparing 
the results of tests on plates and tubes, it may be expected that the tube subjected to f, 
will behave as a flat plate subjected to f,+f,. 


However, sufficient information on flat plates is not available, nor would it be easy 
to compile reliable data for flat plates wide enough for the edge effect to be disregarded. 
Accordingly, it is more profitable to ignore the flat plate and to apply formula (1) to 
compare tubes of several sizes on the assumption that unstable crack length depends 
simply on the attainment of some critical value of f,+/,; or, assuming f, to be kept 
constant, on some critical value of f,/f,. Values of unstable crack length a, computed 
in this way are shown in Table XI for three values of f,//, and for five sizes of tube. 


*See the second paragraph of this Appendix. 


February 1961 31 


H. L. COX AND N. B. OWEN 


The first and third examples relate to tubes tested at Bristol Aircraft Limited, the 
fifth and sixth examples relate to tubes tested at the N.P.L.; the second and fourth 
examples are included to illustrate the variation of a, with 1. 


TABLE XI 


UNSTABLE CRACK LENGTHS (a. IN INCHES) 


f/f 0:5 1-0 2-0 
R=72, t=0-040 12 ie 20 
R=72, t=0-064 18 23 
R=22, t=0-040 5 7 9 
R=22, t=0-064 6 8 10-5 
R=0:-31, t=0-015 0-20 0-27 0:35 
R=0-125, t=0-005 0-075 0-10 0-13 


Appendix II 


RELIEF OF MEAN Hoop TENSION IN A SEGMENT OF SLIT TUBE 


Taking axes Ox circumferentially and Oy axially through the centre of the slit, the 
stress system in the unstiffened tube may be represented in the form 


Circumferential stress, f,= (A cosh zy+Bzy sinh zy) sin 2x. 
Shear stress, q=—X ((A—B) sinh zy + Bzy cosh zy] cos zx. 
Axial stress, f,=— [(A—2B) cosh zy+ Bzy sinh zy] sin zx. 


Here x=n/(2R) and the summations are over all odd values of n. 
The average circumferential stress, f,,, across the axial section at x of the segment 


comprising the slit is 


a/2 x 
> 5 


dy | q dx (A—B) + BC | sin 2x 


where 4=na/(4R), C=cosh 6 and S=sinh 4 


The average value of f,, round the whole perimeter of the tube is 


If we assume f,=0 and f,=f(constant) at y= + (a/2), AC=B(2C—86S) and 
AC + BOS=(4/xn) f, and it follows that 


+1 -T) . Where 7=tanh 6. 
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Evaluation of the series for values of (a/R) from 0 to 1-5 indicates that 


8 R 


Other, perhaps more reasonable, assumptions with regard to the conditions at the 
sections through the ends of the slit would, of course, result in a different relationship; 
but the analysis outlined here may suffice to justify the assumption that the relief of 
mean hoop tension in the way of the slit varies linearly with a/R and to indicate the 
order of the constant of proportionality. 


Appendix III 


Hoop TENSION SET UP AT THE END OF A LONG TUBE BY SHEAR LOADING AND BENDING 
MOMENT APPLIED UNIFORMLY ALL ROUND THE END SECTION 


The differential equation governing the radial displacement y at the section distant x 
from one end of a tube of radius R and thickness t (t < R) under circularly symmetrical 
loading is 


dty 12(l-—o*)y _ 
dx* R?*t? 
of which the solution appropriate to a long tube loaded at one end only is 
y=e-*(A cos 2x+B sin zx), 
where 


Under shear loading F and bending moment M per unit width applied over the 
whole section at x=0, 


A=6(F+Mz2)/(E’2*t?) and 
where E’=E/(1—o”). 
Then, atx=0, y= A 
and the hoop tension f=FEy/R=6(1—o*)(F+Ma2)/(R2&P). 
Or, substituting for z and taking c=0°3, 


f=([2-6F (Rt)-'/7+3-3M (Rt)-') 
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The Vibration of Air in a Duct with a 
Subsonic Mean Flow 


D. S. WHITEHEAD 


(University Engineering Laboratory, Cambridge) 


SUMMARY: The effect of a subsonic mean flow on the organ-pipe type of 
vibration of the air in a duct of uniform area is calculated. The flow at entry 
to the duct is assumed to be at constant stagnation pressure, and a range of 
exit conditions is considered. It is found that the effect of the mean flow is 
to reduce the resonant frequencies by a factor of (1—M7*) (where M is the 
Mach number in the duct) and to introduce damping of the vibration. Except 
when the duct is terminated by an efficient diffuser, this damping is rather 
large and it would require a powerful driving mechanism to produce large 
vibrations in the duct. 


1. Introduction 


The investigation reported in this paper was carried out in connection with a 
suspected case of severe oscillation (not subsequently confirmed) of the air in the 
intake duct of a jet engine installation. It was required to find out what the effect 
of the mean flow would be on the organ-pipe type of vibration of the air in the duct. 
Similar situations occur in the combustion chambers and jet pipes of jet engines, and 
in the combustion chambers of rocket motors. 


The system to be analysed is shown in Fig. 1. A duct of uniform cross-sectional 
area is fed at constant stagnation pressure from an atmosphere through a short 
contraction. The duct is terminated by an area change leading to a large vessel in 
which the static pressure is maintained constant. Small oscillations of the air in the 
duct are considered and, in order to obtain a solution, it is necessary to allow for 
exponential growth or decay of the oscillations. The system is assumed to be one- 
dimensional and fluid friction is neglected throughout. The fluid is assumed to be 
a perfect gas. 


{| 


Ficure 1. Illustration of system analysed. 
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NOTATION 
A,B arbitrary constants 
M-=U/c, mean Mach number in duct 
M. mean Mach number at exit 
T temperature 
U_ mean velocity of flow in duct 
a area of duct 
a, area of exit 
c mean speed of sound in duct 


» Specific heat at constant pressure 


i= /(-—1), indicates component leading 90° in phase 
k constant given by equation (16) 

1 length of duct 

n integer (n> 1) 

p static pressure 

$s a measure of density variation 

t time 


u a measure of velocity variation 


velocity of fluid in duct 
x distance along duct 
2 amplification factor 


y ratio of specific heats of fluid 


logarithmic decrement 
p fluid density in duct 
p, mean fluid density in duct 


» angular frequency of oscillations. 


2 Analysis of Flow in Duct 


The equation of continuity for flow in a duct of constant area 1s 
(pv) rs 0. . (1) 


The equation of motion in the x-direction is 


ov ov 1 op 


Ox pda 
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Since the flow is reversible and adiabatic, and the fluid is a perfect gas, 


and the speed of sound is given by 
YD 
2. (4) 


It is now supposed that p and v vary with time in the following way : 


P=Po {1 see | 


(6) 
v=U {1 +ue“t 


where p, and U are constants, and s and wu are functions of x only and are small 
quantities. The time variation represents oscillations with angular frequency » with 
an exponential growth or decay given by the factor 2. 

Substituting equations (6) into equations (1) and (5), and neglecting terms of 
second order in small quantities, gives 


duds 2+ io 


a+iw @ 


Writing D for d/dx, and eliminating s between these two equations, gives 


(p: u 0 


The solution of this equation is 


(a+hw) 2 (a+bo) 


u=Aei-m © 4 Be (9) 


where M=U/c, and A and B are arbitrary constants. 
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Eliminating the pressure variation from (2) and (3), and using (4), gives 
Aes An 
ct OX Ox 


AIR VIBRATION IN A DUCT 
Substituting (9) into (8) gives the corresponding solution for s, 


AMe 1-™ « +BMe (10) 


The first terms in equations (9) and (10) represent a wave travelling backwards 
at a speed c(1—M), i.e. at sonic speed relative to the fluid. The second terms 
represent a wave travelling forwards at a speed c (1+ M). 


3. End Conditions 


Consider first the end conditions to be applied at the exit from the duct. It is 
supposed that the area change takes place over a very short length of duct, so that 
conditions over this short length can be regarded as though they were steady. 
Conditions at the end of the duct will be denoted by the suffix “2” and at the exit 
after the area change by suffix “e.” Since the static pressure at the exit, p., is fixed 
and the whole flow is isentropic, the state of the fluid at exit is fixed and 7. and p, 
are constants. 


By continuity, 


From energy considerations, 


(12) 


and for isentropic flow of a perfect gas 


3 


Eliminating 7. and T, from equations (11), (12) and (13) gives 
(14) 


This equation is the required relationship between v, and p., since all the other 
quantities appearing in the equation are constants. Substituting for p, and v, from 
equations (6), and picking out terms of first order in small quantities, gives the 
corresponding relationship between u and s in the form 
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If the constant k is expressed in terms of Mach numbers, it is found after some 
manipulation that 


(1+ Me) (14 


(16) 
M.? 


k ‘on 
This is the required end condition at the outlet end of the duct. Three special 
cases are of interest. 

(i) If M.=0 then u=—s/M?*. This is the case when the duct is fitted with 
an ideal diffuser of very large area ratio, so that the stagnation pressure at duct exit 
is constant. 

(ii) If M.=M then s=0. In this case there is no area change at duct exit, so 
that the static pressure at duct exit is constant. 

(iii) If M.=1 then u=4(y—1)s. In this case the duct is terminated by a 
choked nozzle. 

At entry to the duct the assumed condition is one of constant stagnation 
pressure. This corresponds to the special case (i), so that at entry 


. 


4. Solution 


Applying the outlet condition (15) to equations (9) and (10) when x =0 gives 


1+B A+B) 
A | T M 


Applying the outlet condition (15) to equations (9) and (10) when x =1/ gives 


latin _ Gtiv)! (a+iw) ! ore 


Eliminating the ratio A/B from equations (18) and (19) gives 


(1 ~M)(1—kM) 
1—-M'c — 


The real and imaginary parts of this equation must both be satisfied. 
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2al wl (1 —M)(1—kM) 
H 1ce > c(l1—M? 


2al 2wl 


These two equations enable the two unknowns 2 and » to be found. The result 
depends on the sign of the right-hand side of equation (21). 


( (1-M)(1—kM)) 


If | is positive, 
20 
then 27" 
(23) 
2al loo, |. | 
c(l1—M*) (1+M)(1+kM)| 
. 
| 
2wl 
then —M?) (2n—1)= | 
(24) 


In each case n is an integer (n > 1) which determines the various modes of the 
system. 


The solution (23) includes the special cases ((i) and (ii) in Section 3) of constant 
Outlet stagnation pressure and static pressure, and the result is analogous to the 
standard solution for an organ pipe open at both ends. The solution (24) includes 
the special case ((iii) in Section 3) of a choked exit, and is analogous to an organ 
pipe open at one end and closed at the other. 


It will be seen that, in both solutions, one effect of the uniform flow in the duct 
is to drop the frequency from the standard organ pipe value by the factor (1 —M”). 
This simply represents the extra time taken for a wave to travel from one end of the 
duct to the other and back again. 

The other effect of the uniform flow is to introduce damping of the oscillations. 
This is probably best specified by the usual logarithmic decrement 6, which is 
related to the quantity 2 as follows : — 


For the solution (23) this is given by 


(1—M)(1—kM)5S ° 
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For the solution (24) it is given by 


(1+M)(1+kM) 


These equations are illustrated in Fig. 2 (for y=1-4). The full line is the 
boundary between the two solutions (given by KM=1) and also corresponds to 
infinitely great damping of the oscillations. Below this are shown broken lines of 
constant values of né, calculated from equation (26). Similarly, above the full line 
are shown broken lines of constant values of 4(2n—1)6, calculated from equation 
(27). This diagram therefore enables the type of solution to be determined and the 
damping to be estimated for any given values of M and M.. Constant outlet 
stagnation pressure conditions lie on the line OP and the damping is zero. Constant 
outlet static pressure conditions lie on the straight line OQ, and choked outlet 
conditions lie on the line QR. 


R Q 
| 
/| VALUES OF 4 
j A 
= / f 
x 
= | 
Ficure 2. Damping 2 
of duct vibration. = VALUES OF nd r 
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= | oS = 
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DUCT MACH NUMBER M 


It can be seen from Fig. 2 that, except when the exit Mach number is low, the 
damping of the vibrations is high, at least for the lower modes of vibration 
corresponding to low values of n. It is concluded that vibration of this type would 
only occur if a very powerful driving mechanism were present. The exceptions to 
this are cases when the duct is fitted with an efficient diffuser of large area ratio at 
exit, in which cases vibration might be rather easily excited. 
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The Effect of Fuel Distribution on 
Reaction Rates 


M. V. HERBERT 


(National Gas Turbine Establishment) 


SUMMARY: Theoretical reaction rate considerations normally require homo- 
geneity of the mixture supply. An attempt—necessarily crude by reason of the 
analytical difficulties—is made here to cater for the case of non-uniform fuel 
distribution. Assuming a point-source injector in a parallel pipe, the “ loss ” in 
reaction rate is calculated for a limited number of specimen cases. 


1. Introduction 


Since Longwell’s classic experiments with the spherical reactor, much play has 
been made with homogeneous reaction rate theory using effective * global * constants 
derived from his data. Such theories ordinarily assume “ ideal” conditions both 
inside and at entry to the reactor, including the supply of perfectly mixed com- 
bustibles. These conditions could scarcely be satisfied in a practical system, and 
often are not even approached. Heat loss from, and imperfect internal mixing in, 
the reaction zone are inherent in practical designs, and a further departure from the 
ideal can occur in the mixture entry conditions. The fuel may be incompletely 
vaporised, incompletely mixed, or both. The case of non-uniform fuel distribution 
is one for which an attempt may be made to cater analytically in a simple instance. 
Applications of particular interest are reheat systems and ram-jet pilot burners, 
where considerable penalty may be suffered by allowing insufficient distance between 
fuel injection and stabiliser for full diffusion. In such cases, interest is mainly in 
the stoichiometric condition, and it is that to which most attention has been 
given here. 


NOTATION 
A cross-sectional area of duct (ft.*) 
C concentration of fuel in air (Ib./Ib.) 
e turbulent diffusivity (ft.*/sec.) 
P pressure (atmospheres) 
Q air mass flow (Ib./sec.) 


overall fuel flow (Ib. /sec.) 
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radius from centre line of duct (ft.) 


reaction volume (ft.*) 


v ait velocity (ft./sec.) 
vy distance downstream from source (ft.) 
8 weight fraction of air mass flow at given mixture strength 
combustion efficiency 
p air density (Ib. /ft.*) 
equivalence ratio fuel /air ratio 
stoichiometric fuel/air ratio 
Suffixes 
m overall weight mean 
0 centre line 
w wall 
A prime (’) denotes infinite stream conditions. 
Other symbols are defined in the text. 


General Conditions 


The hypothetical arrangement considered is a circular duct containing a single 
central “ point-source ” fuel injector upstream of a flame-holder (Fig. 1). Conditions 
in the duct are assumed to be turbulent and uniform, both as to velocity and 
temperature, and reaction rates are arbitrarily computed for the case of unvitiated 
air at 500°K. In the case of liquid fuel, it is assumed that diffusion and evaporation 
take place simultaneously after injection, and the turbulent diffusivity (¢) is taken to 
be an effective one for the state of the fuel at the station considered. It is reasonable 


i 
i} 
ii “ 
FUEL INJECTOR DIFFUSION PROFILE STABILISER 


Ficure |. Fuel distribution: diagrammatic arrangement of system. 
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to suppose that, as evaporation proceeds, so e will apparently increase with distance 
downstream, until evaporation is complete, when it will assume a constant value. 
In the absence of experimental data for these conditions, we can only take e as fixed, 
i.e. evaporation is assumed to be complete upstream of the flame-holder in 
each case. 


The form of stabiliser shown in Fig. | is not a true representation of either 
a reheat or a ram-jet pilot burner, and we should consider how applicable this 
analysis will be. In computing an overall reaction rate, it is arbitrarily assumed that 
the total reaction volume is shared amongst discrete volumina at different mixture 
strengths so that all reach the same stage of combustion together—i.e. all have the 
same fractional consumption efficiency—and that there is no heat transfer betwee 
them. While this is recognised as a quite artificial treatment, it is the most 
convenient of manipulation in deriving an overall average reaction rate. These 
conditions will give not only a reduction of reaction rate relative to the 
homogeneous case, but also a lower mean temperature. 


These assumptions may be expected to be most nearly realised in a ram-jet 
colander burner, as the stability of the whole reaction there is not dependent on 
propagation from any one particular region and hence on any one set of local 
conditions. It is more difficult to predict how the performance of a reheat gutter 
system will respond to an assumed change in fuel distribution; but qualitatively it 
can be supposed that reduced overall reaction rate will be accompanied by a shift 
in the peak of the stability loop away from the region of stoichiometric overall 
towards some weaker mean mixture strength. 


It is not intended that the quantitative values derived here should be anything 
more than illustrative, bearing in mind the limitations in treatment. The order of 
the effect, however, should be of quite general interest. 


3. Fuel Distribution Profiles 


The radial concentration gradient at any station distant y downstream of a 
point source in an infinite stream is given in Ref. | as 


C O Av [ rev (4ve)] (1) 


QO 4zve 
Other conditions being constant, the concentration at radius r is then 

where k?=4ye/~. 


Suppose a wall concentration C,, be defined at radius ry. Then 


ala 
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A distribution such as that given by equation (3) is indicated in Fig. 1. 


In the practical case of a finite duct, some allowance must be introduced for the 
presence of the duct wall. The method adopted cannot be rigidly defended, but it 
will presently be shown that the result is reasonably satisfactory. The assumption is 
here made that equation (3) can be applied to that portion of the profile inside the 
duct, with r, as the duct radius. Now we know that the total quantity of fuel across 
the section is Q;. If C,, is the overall weight mean concentration, then 


Also, O;= Cpv 2tr dr. 


Hence, substituting from equation (3), 


2(c Ba r dr, 


hich give — (4) 
which gives Cu log (C./C.) 


Now it is convenient to work in terms of equivalence ratio (#) in calculation of 
reaction rates. We know that for the conditions considered 


CO ¢. 
Thus equation (4) may be rewritten as 


Om 
log (o, 
To avoid introducing particular duct dimensions we define arbitrarily, as an 
example, a wall condition of ¢.=0°5 for an overall fuel/air ratio given by ¢,=1 
(stoichiometric). Then equation (5) gives 


1-76 on the centre line. 


Hence, from equation (3), 
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Plotting this relation gives the profile curve (1) in Fig. 2. Other similar curves 
are also shown for values of 9, between 0-5 and 1-0, all for the same overall 
stoichiometric condition. Corresponding to each of these curves in Fig. 2 are those 
of Fig. 3, giving the proportion of air flow by weight (6) reaching any mixture 
strength. 


One must now consider how justified is the assumed treatment of wall effects 
The infinite stream profile (Fig. 1) is modified in practice by collection of fuel on 
the wall and diffusion back into the duct, so that the concentration gradient at the 
wall is zero for equilibrium. For a large duct (or a short distance from the injector) 
this will lead to a distribution as indicated by the upper curve in Fig. 4(a), the thick 
line denoting the infinite stream profile (#’), and the fuel masses represented by areas 
A and B being equal. (Note: On the diagram, area B > area A because of the 
changing radius.) No general diffusion equation can take these wall effects into 
account and, even in a particular case, the treatment of secondary diffusion to give 
the superposed profile (Fig. 4(a)) involves some arbitrary assumptions. The best 
compromise seems to be the retention of a general relation such as equation (1) 
(cf. various empirical curves as, for instance, in Ref. 2) deliberately misapplied 
for convenience. 


WALL 
b= 
4 | dw 
YG J 5 or 4 
CENTRE LINE CENTRE LINE ENTIRE LINE 
(a) True pattern of re-diffu- (b) Substitute assumed for (c) Condition as diffusion 
sion at wall of duct. calculations. nears completion. 


Ficure 4. Effect of duct wall on fuel diffusion. 


The assumption made that equation (3) may be applied to the whole profile 
within the duct gives the upper curve of Fig. 4(b), which may often be a sufficiently 
close representation of the true state of affairs, despite the fact that the concentration 
gradient at the wall is not zero. It would be possible to distort that part of the curve 
adjacent to the wall to give ¢@/ér=0 but, for constant fuel mass flow (OC ¢m) and for 
a given value of @,, this would give a new value of ¢ bearing no relation to that 
computed from equation (5). 


In any case, the final result is only affected if the steps of radius taken in the 
subsequent calculations are sufficiently small and, with our use of 10 per cent 
intervals, it does not seem worth introducing any correction. Concentrating, then, 
on the method of Fig. 4(b) and referring 9’ to the infinite stream profile, the values of 
$w and %,,’ are now compared. 


A 
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For the same overall fuel flow (cx @,,), 


w x 


or dr= | @’rdr. 


Substituting again from equation (3), this yields 


Ow ow) 
(6) 
Applying this to the previous example in which #,=0°5, and for which 


(following Fig. 4(b)) we can write 
= 


we get 4,” = 0°33, which seems a reasonable value. Hence the method is deemed to 
be acceptable for the present purpose and has been used to produce the curves of 
Figs. 2 and 3 already described. 


For a more even distribution, as occurs at a long distance from the injector, it is 
no longer permissible to write 9, = ¢,’, as the picture becomes like Fig. 4(c). This is 
the result of successive re-diffusion from the wall to the centre line and back. In 
this régime there is no real justification for continuing to use part of an infinite 
Stream profile (given by equation (3)) in the manner previously assumed. Since, 
however, the actual shape of the distribution profile becomes relatively unimportant 
as complete diffusion is approached, the method has been retained for convenience 


4. Reaction Rate Calculation 


For a given fuel distribution profile, an assessment of local reaction rate has 
been made over a series of radial steps, using the mass mixture proportions shown in 
Fig. 3. For each step the method of Ref. 3 was used to find the local value of air 
loading. This is defined in terms of the group Q/(VP*), being a function of 
equivalence ratio and inlet temperature, and dependent on combustion efficiency. 
At efficiencies of 100 per cent and close to zero the value of air loading is zero, and 
a peak loading is reached at some efficiency usually around 80 per cent. This peak 
corresponds to combustion extinction or “blow-out”, as stable operation is only 
possible, in general, when the gradient is such that an increase in loading can be met 
by a decrease in efficiency. 


As mentioned earlier, each radial step is treated as a separate adiabatic and 
homogeneous entity without influence on its neighbours, and it has been further 
assumed, purely for convenience, that the same level of combustion efficiency is 
reached simultaneously in each step. Thus a curve of aggregate reaction rate against 
efficiency can be obtained for each fuel distribution pattern. 
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Defining local conditions by 


8=weight fraction of air mass flow at state 
| AV =volume required for reaction of mass 6Q 
Y =air loading of volume Ab Avp:? 


we know that > and = =1-0, 


6Q | Q N 6 
) 
Y 


This result is shown in Fig. 5, corresponding to the distributions in Figs. 2 and 
3, all for ¢,,=1. Fig. 6 gives the blow-out conditions taken from Fig. 5, from which 
it can be seen that in the worst case considered the “loss” of air loading is 
75 per cent. 


Since maximum reaction rates occur close to the stoichiometric condition, and 
fall off rapidly on either side, distributions such as those in Fig. 2 must naturally be 
expected to result in overall losses. It seemed possible, however, that at some other 
value of ¢,, the fractions of mixture nearer the stoichiometric condition might 
balance the effect of those further away. To check this, a single case has been 
worked out for ¢,, = 0°75, with the distribution shown in Fig. 7. The result is given 
in Fig. 8, from which it is clear that the large mass of mixture at conditions less than 
$m has far outweighed the existence of a local stoichiometric region, and that no 
“ gain ” in reaction rate is possible. 


2 4 6C 80 100 20 140 160 


Ficure 5. Reaction rate curves (@,,=1). 
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It is of interest to compare the relative blow-out loading at ¢,,=0°75 (from 
Fig. 8, we have 41 per cent of the uniform distribution case) with an appropriate 
value at %,,=1. Making the assumption that the value of e is unchanged, which is 
experimentally valid for downstream injection, equation (1) leads to 


Hence, scaling the value of 9, from Fig. 7 in the ratio of ¢m, we have that 
at 


at which condition Fig. 2 gives ¢,=0°63. A value of 57 per cent relative loading 
may then be taken from Fig. 6 for comparison with the 41 per cent obtained 
previously. 


The implication is, then, that there is no substantial change in the absolute 
value of “ loss ” throughout the range of overall mixture strength. A stability loop— 
in terms of air loading and equivalence ratio—would be shifted, but not apparently 
blunted. This conclusion is necessarily tentative, because of the crudity of the 
analysis, but may be of some interest. Several practical combustion systems, 
Operating with less perfect uniformity of conditions than Longwell’s spherical 
reactor, have produced stability loops much blunter than this, and differences in 
fuel distribution might have been a contributory cause. 


5. Conclusion 


A tentative method has been used to plot the fuel distribution in a finite duct. 
The reaction rate has been calculated for a number of cases with an overall 
stoichiometric mixture, and the “loss” due to uneven distribution is found to be 
75 per cent in the worst case considered. From limited work at other conditions, it 
is thought that no significant change in the shape of a stability loop should result 
from fuel distribution alone. 
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The Behaviour of a Clamped Circular Plate 
in Compression 


A. N. SHERBOURNE 


(University Engineering Laboratory, Cambridge) 


Summary: A theoretical solution is presented for the problem of the clamped 

circular plate loaded in uniform compression. The solution employs a numerical 

method programmed for a digital computer. Instead of solving the classical 

von Karman large deflection equations, a step-by-step integration of the elastic 

differential equations of equilibrium is carried out until suitable boundary 

conditions are attained. The method is an extension of one developed earlier 
to explain the behaviour of the simply-supported plate. 


1. Introduction 


The ultimate load for thin plates under uniform compression may be estimated 
from the intersection of an elastic-plastic loading path and a rigid-plastic mechanism 
line, when the resulting load will represent an upper bound on the collapse of the 
plate. The treatment of the circular plate appears to fall naturally into two distinct 
problems, involving the simply-supported and the clamped plate; each problem may 
be regarded separately by defining the complete loading history for plates with these 
particular boundary conditions. For the plate with the elastically restrained edge, 
an estimate of the collapse load and loading history may conveniently be obtained 
by interpolating between the two solutions mentioned. 


The classical problem of the clamped circular plate in compression has been 
discussed by Bodner“. His solution applied power series and perturbation methods 
developed earlier by Friedrichs and Stoker to explain the behaviour of the simply- 
supported circular plate loaded beyond the critical thrust. Timoshenko“ has used 
numerical methods to determine the elastic behaviour, at large deflection, of a 
circular plate bent by couples uniformly distributed along the supported edge. 


A previous paper’ developed a numerical procedure for the post-buckling 
behaviour of a simply-supported circular plate. The method obviated the solution 
of the classical von Karman large deflection equations, substituting, instead, a step- 
by-step integration of the equations of equilibrium until specific boundary conditions 
were attained. Such a method had the added advantage of accommodating plastic 
behaviour by a suitable modification of the equations of moment, thrust and 
equilibrium for the bent plate. 
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The current investigation applies this numerical method to the clamped circular 
plate loaded in uniform compression beyond the critical thrust, the analysis being 
confined to elastic behaviour throughout the plate. The procedure converts the 
conventional two-point boundary-value problem, in which statical and geometrical 
conditions are fulfilled at the centre of the plate and along the supported edge, into 
an initial-value problem by selecting appropriate numerical values for the moment 
and thrust at the centre of the plate. These values completely define the geometry 
of the plate as expressed in terms of the radial and transverse displacements and are 
so chosen as to preserve elastic behaviour. A step-by-step integration of the 
differential equations of equilibrium is carried out and statical and geometrical 
conditions are evaluated over each interval. This process is continued until a 
critical radius is reached when the statics and geometry match prescribed boundary 
conditions on the problem. For the clamped plate this implies a condition of zero 
slope, and the boundary tractions are evaluated in this position. By varying the 
initial combinations of moment and thrust, a set of boundary tractions can be 
obtained corresponding to plates of various size. For a plate of specific radius, with 
particular boundary conditions, the elastic, post-buckling behaviour at large 
deflection may be obtained by interpolating between the contours of moment 
and thrust. 


NOTATION 
u,w_ radial and transverse displacements 
r radius 
R_ plate radius 
Poisson’s ratio 
h_ plate thickness 
D flexural rigidity of plate, = Eh’ / {12 (1 — 
E modulus of elasticity 
y.8 non-dimensional moment and thrust 
y,,8, non-dimensional moment and thrust in radial direction 
y:.,8, non-dimensional moment and thrust in tangential direction 
m, fully plastic moment per unit width of plate, = h*f,/4 
fy yield stress in tension or compression 
o mean stress in plate 
kK? = E/ {3f, (1 
a=dw/dr 
q’=hp/D 
p applied edge load per unit width 
P.r elastic critical load per unit width 
y, B=y/k* and 8/k* respectively 


w,u,r, R=w/h,u/h,r/hand R/h, respectively 
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W,U_ functions of 
P total edge load 
P, yield load. 
The suffixes r, t, 0 denote the radial and tangential directions and the centre of 
the plate respectively. 
2. General Analysis 


The elastic critical load and subsequent post-buckling behaviour of a circular 
plate may be obtained from the equilibrium of the plate at large deflection, 


dr’ 2r ) dr )( dr’ ) 
1 1 dw 12 du (dw 


where u and w are the radial and transverse displacements respectively at radius r, 
u 18 Poisson’s ratio and h is the plate thickness. 


The critical load is obtained from equation (2) by considering the associated 
plane strain problem in the small deflection theory 


da dz 
3 

such that z=dw/dr, g’?’=hp/D, where D is the flexural rigidity of the plate, 
Eh* / {12 (1 — »?)}, and p is the applied edge load. 


For the clamped plate, of radius R, a solution of equation (3) is given by the 
characteristic value 


14-68D 
Pes 


The moments and thrusts obtaining in the plate can be expressed, in terms of the 
parameters y and 8, as 


yumm, 


B=a/f,, 


where m is a general moment in the plate, m, the full plastic moment per unit width 
of plate (=4h°f,), « the mean stress in the plate and f, the yield stress in simple 
tension or compression. At large deflection, the moments and thrusts are determined 
from the system of strains and stresses obtaining in the bent plate“ and may 
be written as follows. 
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k2} E dw 
u du dw\?* ] 
[+ ar +5 (4) 
where k?= 3f, (1 


and the suffixes r and ¢ denote the radial and tangential directions respectively. 


The initial conditions are obtained from the symmetry of statics and geometry 
at r=0. Substitution for identical curvatures, strains and slopes in the radial and 


tangential directions in equations (1), (2) and (5), gives 
dw 
(6) 
Ldr 3k? (1 +p) 
| =0 
dr’ 
and y,=7,=y7,. 
A limiting condition of elasticity throughout the plate is derived by 


considering the interaction between moment and thrust within a maximum stress 


criterion 


such that 


(7) 


3, caters for membrane tension or compression. 


The previous equations may be non-dimensionalised through the relationships 


when the entire problem can be defined as follows. 
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(a) The equations of moment and thrust, 


dr’ 


Y=" 


dw 
r dr 


a(S) + | 
dr 
(8) 
r dr dr 
du. (dw? 
(b) The differential equations of equilibrium, 
_ au du u_ (dw dw 
"dP dr? 2 * 
9 
; dw dw 1 dw. 127 1 / ©) 
(c) The initial conditions, 
a + pt) 
(10) 
du 
i 
& 7 
dr 
(d) The conditions for elasticity. 
Elastic behaviour at the centre of the plate can be expressed by the inequality 
a 


whereas at radius r, the inequality 
directions is given by 
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3. The Numerical Method for Integrating the Differential Equations 
of Equilibrium 


Before suggesting a programme for the numerical integration of the differential 
equations of equilibrium, it is necessary to obtain finite initial values for d*u/dr’ 
and d°w/dr’. The singularities in these functions, which prevail at the origin, are 
of the form 0/0 and may be removed by introducing the transformations 


w (r)= W (r’) 


when the entire problem may be redefined in terms of the new variables, 
as follows. 


(a) The equations of moment and thrust, 


y-= W'+4xW"] 
B,= 3[(1+ ») U+2x | 
yer (14) 
—(2(1+ 4uxW"] | 
B.= {U’+(W}] 
where primes (’) denote differentiation with respect to x. 
(b) The equations of equilibrium, 
u—3 
2 (wy - | 
x 2x 
; (15) 
= — W" + + OW" (WY) 
(c) The initial conditions, 
B 
W J 
U 3 ( 1+ 
WwW,’ Yo — 3)y 
2(1+ 4) Uo 16 
16 
Ww,” Y Fe U 5) Boy." 
4(14+ 24 (1+ p)? 
0 
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(d) The yield conditions, 


~ 
iN 
w 
~ 
t 


(18) 


The load-deflection relationship for the clamped circular plate in compression 
can now be obtained using a numerical procedure, the equations (14) to (18) being 
of a form suitable for programming on a digital computer. The various stages of 
such a programme may be outlined as follows. 


Stage 1. Select values of y,, ©, in accordance with equation (17) 

Stage 2. Evaluate the initial conditions as given by equation (16) 

Stage 3. Advance the functions W, W’, W”, W’’, U, U’, U”, a single step 6, using a 
numerical integration procedure, and evaluate the moments and thrusts as given by 


equation (14). Check the yield conditions, equation (18), to ensure that the plate is 
elastic at radius 6. 


Stage 4. Repeat Stage 3 until the slope is zero. This implies that W’=0. 


Stage 5. Evaluate and print out the required results. For the clamped plate, 
obtain 


the critical radius R, the boundary moment yz, 
the boundary thrust {,=p., the total deflection w,,. 
‘tage 6. Repeat Stages | to 5 using a different set of y,. 8). 


4. The EDSAC Sub-Routine 


The equations (14) to (18) were programmed for EDSAC II, the Cambridge 
University High Speed Digital Computer, the numerical integration being carried out 
using a modified fourth-order Runge-Kutta process which is incorporated in the 
machine as a Library Sub-Routine. As contained in EDSAC, the sub-routine will 
handle a set of simultaneous first-order differential equations in which each 
derivative is expressed explicitly in terms of the dependent variables 


Vn) 


February 1961 57 


aft 
B, )| 


A. N. SHERBOURNE 


Any equation or set of equations must be expressed in this form before the 
integration process can be applied. 


Since the differential equations of equilibrium are of the form 


= F (x, Ww’, w”) 
= G (x, WwW’, WwW”), 
they can be reduced to a set of first-order equations by introducing as many 
dependent variables as may be required, 
X=5, 
U’ =s, 
Ww =, 


ds, 
Ww’ = 

The resulting set of simultaneous equations is now of the form required in the 


sub-routine 


ds, 
F (s,, 535 S,) 
ds, 
=G U, Ss» S4) 
dx 
dx 
ds, 
dU 
dx 


The sub-routine operates on the equations (19) by evaluating the dependent 
variables s,, U, s,. s, and s, at the half-interval 6/2 using a fourth-order Runge- 
Kutta process and uses these functional values at the half-interval to advance the 
functions over a full interval. In this respect the procedure differs from orthodox 
numerical integration techniques where functional values are obtained directly over 
the full interval 6. 


The Runge-Kutta procedure, used in the sub-routine of EDSAC for the 
numerical integration of a set of first-order simultaneous differential equations, has 
powerful advantages when compared with other methods of numerical integration 
inasmuch as it can easily be modified to include changes in the step interval, it is 
self-starting and does not require back information to continue the integration 
procedure. These features become important in high-speed digital computers since 
they lead to a minimum of special instructions, with consequent saving in 
machine time. 
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APPLIED EDGE LOAD f= p/k 
FiGuRE 1. Plate size against applied edge load for a clamped circular plate (for elastic 
behaviour in the plate). 
Notes: B Xk? Xhf, =mean stress at the centre of the plate (tons/in.) 


ot al * =moment at the centre of the plate (ton in./in.) 


p*k* X hf, =edge load (tons/in.) 


k?= 322 for E= 13,400 tons/in.*, f. =15:25 tons/in.2, «=03 
5. Results 
The programme was evaluated for the range of 7, and 2, given in Table I, 


positive 8, being associated with a tension field and negative 8, with a compression 
field at the centre of the plate. The relationships between plate size, edge thrust, 
edge moment and central deflection are given in Figs. 1, 2, and 3 through contours 
of the moment and thrust parameters acting at the centre of the plate. For a given 
plate size, the load-deflection relationship may be determined by interpolating 
between contours. In deriving the results, the values chosen for the mechanical 
properties of mild steel were 


E= 13,400 tons/in.’, 
f,=15-25 tons/in:’, 


k? = 322, and hence 1/k?=3-1 x 
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Ficure 2. Plate size against edge moment for a clamped circular plate (for elastic behaviour 
in the plate). 


6. Conclusions 


A general load-deflection relationship for the clamped circular plate in 
compression may be deduced from the contours of Figs. 1 and 3. It can be shown 
that the relationship is of the form (see the Appendix) 


where the applied edge load and central deflection are expressed non-dimensionally. 
It should be noted that the first term of the expression on the right hand side 
represents the critical load as given by equation (4) and the second term represents 
the elastic, post-buckling behaviour. Equation (20) may be fitted to the results 
obtained from the computer programme by plotting the load against the square of 
the deflection for plates of various sizes when the locus of all points, for a given 
plate size, will constitute a straight line intersecting the load axis in the critical load 
and having a slope which defines the coefficient —k,(E/f,4h/R)*. Using this 
approach, the coefficients of equation (20) may be approximated to 


k,=1:348 and k,=0-22. 
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Figure 3. Plate size against central deflection for a clamped circular plate (for elastic 
behaviour in the plate). 
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Appendix 


The form of th. load-deflection curve near the elastic critical load may be deter- 
mined from the totai energy of the bent plate, 


V=V, +. Vu +V;, (21) 
where the total energy V is made up of a bending component V,, a membrane com- 
ponent V,, and the work done by the external load Vy. For large deflections the total 
energy can be written in terms of the total applied load P and the radial and trans- 
verse displacements u and w respectively. 

V=¢, (w*)+9, (u’, uw’, w') +9, (Pu), (22) 
such that V (w’), 
Vy =>. uw*, w*), 


V (Pu). 


From the theorem of minimum potential energy 


eV 
=() 
cw 
and u may be obtained as a function of w, 
The total energy is now a function of w only, 


The load-displacement relationship may be obtained from the condition 


=0=V, (w)+v, (w*)+¥, (Pw), : (25) 


when P=P..,.+ Kw’. 
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The Design of Wing Plan Forms for 
Transonic Speeds 


R. C. LOCK, M.A., Ph.D. 


(Aerodynamics Division, National Physical Laboratory) 


SUMMARY: A method is suggested for designing the shape of the plan form of 
an unwarped swept-back wing so as to control the magnitude of the singularity 
that occurs in the chordwise loading distribution at the leading edge. At a 
Mach number of unity it is possible to do this simply and directly by an 
approximate application of linearised lifting-surface theory. Details are given 
of three families of wing plan forms, of varying sweep and aspect ratio, in all of 
which the strength of the leading-edge singularity is held constant near the wing 
tips, behind the Mach line from the root trailing edge. The local chordwise 
load distributions have been calculated in detail for several of these wings and 
it is found that in all cases the loading curves become effectively independent 
of spanwise position over that part of the span for which the singularity 
is constant. 


1. Introduction 


In the design of swept-winged aircraft to fly at transonic speeds, one of the 
important features which must be chosen is the wing plan form. The aspect ratio 
and basic angle of sweep will probably be fixed within narrow limits as a 
compromise between aerodynamic and structural considerations, but the detailed 
shaping of the plan form can still have a considerable effect on the aerodynamic 
behaviour of the wings. This is true to some extent for symmetrical wings at zero 
incidence, but the lifting aspects of the problem are much more important and are 
all that will be considered here. 


The treatment must be restricted at the present time to the use of linearised 
theory for lifting surfaces, so that the effects on the loading of the wing thickness 
distribution, as well as of viscosity and shock waves, have to be neglected. Even so, 
at supersonic speeds there does not seem to be any satisfactory method applicable to 
wings of arbitrary plan form, though a great many particular solutions are available 
for wings with straight leading and trailing edges‘, and the methods which have 
been developed for subsonic speeds‘ all break down as the Mach number 
approaches unity. Fortunately, however, at M=1, lifting-surface theory is con- 
siderably simplified and becomes identical to “sender wing” theory®’. This form 
of the theory will be used in the present paper. It is felt that this should not unduly 
restrict the practical applications, since the results for highly swept wings should be 
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directly relevant over a range of Mach numbers around unity, and even outside this 
range the conclusions should be qualitatively correct. 


There is now some experimental evidence that this theory does give good 
results for overall wing lift, provided that the leading-edge sweepback exceeds about 
45°. Further, a comparison with theory of the detailed loading obtained at low 
incidences in a recent series of tests at the N.P.L.” on an uncambered wing of aspect 
ratio 2:8, taper ratio 1/3 and half-chord sweep 45° has shown reasonable agreement, 
particularly over the inboard half of the wing. This comparison is shown in Fig. 1. 
It can be seen that the agreement is poorer near the tip; this is partly due to the 
complicated shock system which appears in this region at quite small incidences 
(see Fig. 1) and to the failure of the theory to predict correctly the position of these 
shocks. However, this shock system, and particularly the “forward shock”, which 
does not occur at all in the theory, is to some extent associated with the very high 
suction peaks and local Mach numbers which occur near the tips of conventional 
tapered swept wings at incidence in the transonic speed range. These high suction 
peaks are a fundamental defect of wing plan forms with straight edges, particularly 
when highly tapered, and it was felt that if they could be alleviated by suitable plan 
form design—and this is one of the principal objects of the present paper—then not 
only would such plan forms be of some practical value, but might also show better 
agreement between experiment and theory. 


The present method of design is similar in some ways to that used by Brebner“? 
in designing “constant C,,” plan forms for low speeds, but concentrates more directly 
on the real source of the trouble which, in linearised lifting-surface theory, appears 
as the singularity* in the load distribution at the leading edge. The trailing edge is 
first chosen to be straight and uniformly swept (see Fig. 2); this is not an essential 
feature but it represents the most feasible shape in practice if conventional controls 
are to be used. At the centre section of any swept wing (or at the wing-body 
junction of a complete aircraft) the load distribution has no singularity, so that over 
the inboard part of the wing the strength of the leading-edge singularity inevitably 
increases steadily in the spanwise direction; the rate of increase depends on the 
shape of the inner part of the leading edge, which may be chosen in various ways, as 
described in Section 2.4. The essential feature of the present method is that, out- 
board of a certain spanwise position, which may also be chosen arbitrarily, the 
leading-edge shape is modified in such a way as to keep the strength of the 
singularity constant out to the tip. The position P, (see Fig. 2) at which the leading 
edge modification begins determines the final aspect ratio of the wing. Thus, for a 
given trailing edge sweep and initial leading edge shape, a single family of wings of 
varying aspect ratio may be developed by varying this position, or alternatively by 
keeping P,, fixed (in the streamwise direction) and varying the initial leading-edge 
shape. Three such families are considered in the present paper, each of which is 


*The strength of the leading-edge singularity is one of the factors which influence the magnitude 
of the suction peak in real flow, although this must also depend on other factors, notably the 
thickness distribution of the wing section near the leading edge. 
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FiGuRrE 1. Comparison between sonic loading theory and experiment: M=1-0. 
Aspect ratio = 2-828 Taper ratio=1/3 Wing section R.A.&. 102, 6% thick 
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believed to be of some practical interest. In all cases the method leads to a 
characteristic curved tip shape, as shown in Fig. 5 (p. 77). 


An interesting feature that is apparent when calculations of the detailed loading 
are made for these wings is that, once the strength of the leading edge singularity is 
held constant, the shape of the chordwise loading curves also becomes almost 
independent of spanwise position. It follows that if the thickness distribution of 
such a wing is designed to give fully swept isobars at zero incidence, then this 
property should be retained in the tip regions under lifting conditions; and similarly 
for a wing, cambered and twisted to give fully swept isobars at one particular 
lift coefficient, under off-design conditions. 


NOTATION 
x,y,z Cartesian co-ordinates, using root chord as unit of length 

c (wing chord)/(root chord) 

s (semi-span)/(root chord) 

S (wing area)/(root chord)’ 

aspect ratio 

M Mach number 

a angle of incidence 


% angle of sweepback of leading edge 
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¢, angle of sweepback of trailing edge 
€=(x-—x,)/c 
n=ytan¢ 
C, pressure coefficient 
1 loading coefficient =C, (lower surface) —C, (upper surface) 
k= y, (x)/y: (x) 
/(1—k?’) 
x= E (k’)/ K (k’) 
K (k’) complete elliptic integral of the first kind with modulus k’ 
E(k’) complete elliptic integral of the second kind with modulus k’ 


| (82%) 


H(x) function determined by equations (5) or (6) 
o strength of leading-edge singularity (see equation (7)) 
7, constant value of « in tip region 
m=(dx,/dn),= (leading edge slope)/(trailing edge slope) 
G =H? (1 —k*) (equation (15)) 
A parameter used in Section 2.4.3. 
C.(y) local section lift coefficient 
C, wing lift coefficient 


§ Spanwise loading factor 
i 
K_ vortex drag factor (equation (28)) 
Cy pitching moment coefficient about x= 0, based on root chord 


xX position of centre of pressure relative to root chord 
Suffixes 

| leading edge 

t trailing edge 

1 value at x=1 


0 value at P, (Fig. 2). 


2. Details of the Design Method 


SUMMARY OF SONIC LIFTING-SURFACE [THEORY 
The co-ordinate system and basic notation used are explained in Fig. 2. Axes 
Ox and Oy are taken in the usual way, with the origin O at the apex of the wing, Ox 


February 1961 69 


in the streamwise direction and Oy in the spanwise direction. It is convenient to 
take as unit of length the root chord c,, and to use in the spanwise direction the 
co-ordinate 


n=y tang, ‘ (1) 


where ¢, is the constant angle of sweep of the trailing edge; thus in the (x, n) system 
of co-ordinates the trailing edge has the equation 


If y, and y, are the y co-ordinates of the leading and trailing edges respectively 


at a given value of x, then the loading coefficient /= AC, for a flat wing at incidence z 
at M =| is given by 


l yi dy, | 
4x /(y,? — dx 
or ™ dn 
ltang, dx 
and by 
l tan 9, dy 
H (x) cot /(1 —k” sin? x) 
X 
+ (1 —k’* sin? y)d -«| dx for x>1 (4) 
y, (x) (x) 
h = : 
Q) 
k’= /(1 —k?), 
_ E(k’) (K and E being complete elliptic integrals of the first 
K (k’) and second kinds respectively), 
and =sin ="). 


The function H (x) has to be determined from the fundamental integral equation 


dn, =0, 5 


where »), is the value of », at x=1, and », is considered as a known function of », (see 
Ref. 2, equation (50)). 
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The following alternative form (Ref. 2, equation I1.5) of equation (5) has proved 
useful for checking values of H obtained approximately by other methods: 


dn. 


_ {H(n)—H 


(n? 2 mv dn). (6) 


Here H (1) means the value of H at the point at which the leading edge has the 
co-ordinates (x, 7). 


2.2, THE LEADING-EDGE SINGULARITY 


Near the leading edge the elliptic integrals (the second and third terms) in 
equation (4) are negligible, and it may easily be shown that the loading in the 
chordwise direction has the singularity given by 


tan 


Ve) as £—>0, (7) 
a 
where €- 
is the fraction of the local chord measured from the leading edge x,{), and o is 
given by 
” 
if 
(8) 
it 
c dx, 
x, (y)-1 
Note that k= = (9) 
and c(n)=1+—x, : (10) 


Before proceeding further, it is of interest to calculate the strength « of the 
leading-edge singularity in some simple cases. 


(a) For an infinite swept wing of constant chord, 


d 
for which gus =I, 
dx; 
we have for nx | 


+ (11) 
Vn for »>1] 


February 1961 71 


R. C. LOCK 


ASYMPTOTIC VALUE=V2 
=. 
Tre) 


2 3 4 5 
” 


FiGuRE 3. Leading edge singularity for an untapered swept wing. 


where H has been calculated by Mangler” and Mirels®. These values of « are 
shown plotted against 7 in Fig. 3. It can be seen that o increases steadily in the 
spanwise direction, but approaches the asymptotic value /2 (corresponding to an 
“infinite yawed wing”) quite slowly; for example, when »=5 (aspect ratio = 10 cot 4) 
* has only reached 90 per cent of the final value. 


(b) For a highly tapered “swallow-tail” wing with leading edge x, = #7, it is 
known®: » that H/(1 —k*) ~ 1 when m > 1, 


so that Cm for all values of ». 
vim Cc 
Here c=i-(m-1)», 
vm 
(1 + § tan 1 —(y/s) 
| semi-span 
where s= is the rati 
(m—1)tang, chord 
V(1+stan®@,) . 
The quantity o “is plotted against y/s in Fig. 4. 


s tan 9, 
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FIGURE 4. Leading edge singularity (~) for a highly tapered wing 


The most noteworthy feature in Fig. 4 is the rapid increase in o as the tip is 
approached. This is clearly an undesirable feature of such wings, and invariably 
leads in practice to flow separation at quite small incidences. For cropped wings of 
this type with straight streamwise tips, the leading-edge singularity takes exactly the 
same values as at the corresponding points on an equivalent fully tapered wing, so 
that, for wings with taper ratios greater than about 4, the increase in © is less 
serious, but there is then a loss of lift in the tip regions* and additional shock waves 
may appear. In the next section it will be shown how, by designing curved leading- 
edge shapes so that « remains constant towards the tips. both these disadvantages 
may, in principle, be avoided. 


2.3. Trp SHAPES WITH CONSTANT LEADING-EDGE SINGULARITY 


It is clear from the theoretical treatment summarised in Sections 2.1 and 2.2 
(and also from physical considerations) that at M =| the strength o of the singularity 
at a particular point on the leading edge depends only on the shape of the leading 
edge upstream of that point. It is, therefore, possible in principle to specify 7 as 
any reasonable function of x (or 1), and to design the leading edge step by step so as 
to give the desired values of «. In all the cases that will be considered here we 


*Linearised theory for M=1 predicts zero loading in the region behind the Mach line from 
the tip leading edge'?); such regions of zero (or even negative) loading are also observed 
experimentally'*), although their boundaries are not predicted correctly by theory 
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shall take > to be a constant o, whenever x => 1 (downstream of the Mach line from 
the root trailing edge); this assumption appears to give sensible and useful tip 
shapes. Upstream of this line (x < 1) there are many alternative possibilities, some 
of which will be considered in Section 2.4. 


When o=c, is constant for x > 1, equation (8) can be rewritten to give the 
following equation for the leading-edge shape 


ay _ (1+ 
dx (1—k?) 
‘ 2 
H?(1+k)’ 
here the suffix / has been dropped and k=(x—1)/n. It is convenient to use k as 
independent variable; then 


dyn _ | _ (x—1) dk 


dx k ds 


so that equation (13) becomes 


dk o,? k? 
(x-1) =k- (14) 
If we write H? (1-k?)=G : ai GES) 


(the reason for this will appear later), equation (15) can be arranged in the form 


dx =ak| ] 
kK 


which can be integrated to give 


t 
where », is the value of » when x=1. 


If G (or H) were known as a function of k, equation (16) would provide an 
explicit solution to the problem. In fact, H is not, of course, known explicitly, but 
depends on the complicated integral equation (5) (or (6). The procedure adopted in 
the present paper is to obtain a suitable approximation for H and to use it in 
equation (16) to find a first approximation to the leading-edge shape. A better 
approximation to H can then be found by means of equation (6), and this could be 
used in (16) to improve the leading-edge shape; but in practice this has not proved 
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necessary (see Section 2.5). A suitable approximation to H is provided by Mirels 
(Ref. 3, Appendix D). He shows that 


where m= dx/ dn, is the ratio of the trailing-edge slope to the leading-edge slope for a 
given value of x. The form of the terms neglected on the right-hand side indicates 
that equation (17) should provide an adequate approximation to H provided either 
that k is small or m is large. In fact Mirels®’ has shown that, for an untapered wing 
(m=1), the approximation is valid as far as k=4; while Mangler’s results indicate 
that equation (17) may be used for all values of k provided that m is greater than 
about 2. Now all the plan forms developed here have the property that m is quite 
large wherever k is greater than 4, and indeed becomes infinite at the tips (k = 1), so 
that it seems reasonable to use equation (17) throughout. This enables an explicit 
approximation to H to be found, as follows. 


Using equation (13), 


dx H’(1+k) 


we have 
dy 


Substituting this expression for mm in equation (17), 


dk 
we get ) +k) log (4/k) (18) 
Writing HJ (1—k*)=z(k) (so that G=z’*). 


and differentiating equation (18) with respect to k, 


dk log (4/k) 


we get 2? 
which may be integrated to give 


(1 —t)dt 
30.2 on k=0) 
log (4/1) (since 1 when 


{ Ei 4 Ei (210g) } 


x 


lt 
where Ei (log x)= 
og 
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TABLE I 
VALUES OF THE FUNCTION G (k, 7 ,) 
2=0'5 0-6 0-7 0-8 0-9 1-0 
: k=0 1-0 1-0 1-0 1-0 1-0 1-0 i 
0-1 0-9790 0-9748 0-9705 0-9663 0:9621 0-9578 | 
0:2 0-9528 0-9432 0-9336 0-9239 0-9141 09044 | 
0:3 0-9249 0-9096 0-8940 0-8784 0-8626 08467 
0:4 0-8970 0-8751 0-8542 0-8323 0-8102 0-7877 
0-5 0-8703 0-8432 0-8156 0-7876 0-7591 0-7300 
0°6 0-8458 0-8133 0-7802 0:7463 0-7116 0-6760 
0-7 0-8248 0-7875 0-7494 0-7102 0-6699 0-6284 
0:8 0- 8083 0-7673 0-7257 0-6816 0-6368 0-5902 
0-9 0-7976 0-7540 0-7091 0- 6628 0-6148 05649 
TABLE II 
LEADING EDGE SHAPES: VALUES OF y/y 
7, 2=05 0-6 0-7 0-8 0-9 1-0 
k=0 1:0 1-0 1-0 1-0 1-0 1:0 
0-1 1-0503 1-0612 1-0722 1 -0836 1 -0952 1-107] 
0-2 1-1010 1-1242 1- 1487 1-1745 1-2017 1 +2305 
0-3 1-1509 1-1880 1-2283 1-2720 1-3197 1 -3720 
0-4 1- 1987 1 +2509 1 -3089 1-3740 1-4476 | 
0-5 1-2428 1-3102 1 -3872 1-4762 1 -5807 1-7053 
0-6 1-2815 1 - 3633 1-4589 1-5725 1-7104 1-8820 
0-7 1-3129 1-407] 1-519] 1-6552 1-8248 2-0434 
0-8 1-3358 1 -4392 1 -5636 1-7170 1-9116 2° 1683 
0-9 1 +3495 1-4583 1-5901 1-9631 2°2424 
1:0 1-3539 1-4644 1-5985 17652 1-9791 2+2650 


Thus we obtain finally the approximation* 


G =H? [1+ Bi (log 4) —48i (210g 4)} ] 


The functions G (k, o,) given by equation (19) are tabulated in Table I for «,?=0'5, 
0:6, 0:7, 0:8, 0-9 and 1:0. The corresponding values of 


(1 ~k) dk 
=exp | 


have been found by numerical integration, and are given in Table II. 


When », and o,? are fixed, the shape of the leading edge downstream of x= | 
can then be found parametrically in terms of k from equation (16) together with 


x=1+kn (equation (9)). 


*Note that this is only valid for those tip shapes for which « is constant. 
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The curves thus obtained depend only on the two parameters 7, (the value of » at 
x=1) and o,? (the value of dn/dx at x=1); but for a given value of o, the curves 
obtained by varying 7, are geometrically similar with respect to the point (1, 0). 


2.4. THE FORWARD PART OF THE LEADING EDGE 


In order to fix the plan form precisely, it is now necessary to specify the shape 
of the leading edge upstream of the Mach line (x=1) from the root trailing edge. 
This may of course be done in an infinite number of ways; in the present section we 
shall consider only three families of shapes which appear physically sensible, but 
other such families could easily be constructed. 


2.4.1. Type A 


The simplest course is to take the leading edge to be straight over the whole 
segment OB (Fig. 2). Its equation will then be 


x=m,n (x1), 


where m, is the ratio of the slope of the straight portion of the leading edge to that 
of the trailing edge. In this case we have clearly 


and (from equation (8) or (13)); 


0 


with these values of 7, and «,* the tip shape may then be found as described at the 
end of Section 2.3. Each value of m, thus defines a unique plan form (apart from 
the trailing edge sweepback ,, which may be chosen arbitrarily). Three typical 
examples corresponding to 7,*=1/m,=0-5, 0-7 and 1-0, are shown in Fig. 5(A); the 
ordinates of the leading edge are given in Table III (Type A). In each case ¢, has 


TABLE Ill 
COMPLETE LEADING EDGE SHAPES 
Type A: Leading edge straight forx <= 1 (n=0,?x) 


k x x x 

0 1:0 0:5 1-0 0:7 1:0 1-0 
0:1 1-0525 0-5252 1-0751 0-7505 1-1107 1-1071 
0:2 1-110] 0-5505 1-1608 0:8041 1-2461 1-2305 
0:3 1-1726 0-5754 1-2579 0-8598 1-4116 1 -3720 
0-4 1 +2397 0-5994 1 -3665 0-9162 1-6127 15317 
0:5 1-3107 0:6214 1-4855 0:9711 1-8526 1-7053 
0-6 1 -3844 0-6407 1:6127 1-0212 2-1292 1 -8820 
0-7 1-4595 0-6565 1-7444 1-0634 1-4304 2-0434 
0-8 1 -5343 0:6679 1:8756 1-0945 2-7347 2-1683 
0:9 1 -6073 0:6748 2-0018 1-1131 3-O181 2°2424 
1-0 1-6770 0-6770 2-1190 1-1190 3-2650 2:2650 
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TABLE IlI—continued 
Type B: Leading edge straight for x <= (n=x) 
(The plan form for 7 7=1-0 is identical with the 
corresponding plan form of Type A) 


Os 2—()-7 
x 
—1-0 0-5 0-5 
0- 5266 0-5256 
0-5569 0-5539 
—0-7 0-5912 0-5840 
0-6303 0-6162 
—0-5 0-6747 0-6506 
0-7252 0-6871 0-7143 0-7142 
0-7823 0-7257 0-7698 0-7675 
—Q-2 0-8468 0-7659 0-8353 0-8238 
—0-1 0-9192 0-8076 0-9115 0-8851 
0 1-0 0-8499 1-0 0-9501 
0-1 1 -0893 0-8927 1-1019 1-0187 
0-2 1-1871 0-9357 1-2183 1-0913 
0-3 1 -2934 0-9781 1-3501 1-1669 
0-4 1-4075 1-0188 1-4974 1 -2435 
0-5 1 -5282 1-0563 1-6590 1-3179 
0:6 1-6535 1-0892 1-8316 1-3860 
0:7 1-7811 1-1159 2-0103 1-4432 
0:8 1-9083 1-1354 2:1818 1-4855 
0-9 2-0323 1-1470 2-3596 1-5107 
1-0 2-1507 1- 1507 2:5187 1 -5187 
Type C: Leading edge curved throughout 
(The plan form for ~ 2=1-0 is identical with the 
corresponding plan form of Type A) 
Os 0 
k x x 
0 0 0 0 
0-1 0-0759 0-1 0-0869 
0-2 0-1499 0-2 0-1727 
0-3 0-2217 0-3 0-2571 
0-4 0-2914 0-4 0-3401 
0-5 0-3586 0-5 0-4215 
0-6 0-4232 0-6 0-5013 
0-7 0-4850 0:7 0-5792 
0-8 0-5437 0-8 0-6551 
0-9 0-599] 0-9 0-7288 
0 1-0 0-6510 1-0 0-8001 
0-1 1 -0684 0-6838 1-0858 0-8578 
0-2 11-1433 0-7167 1-1838 0-9190 
0-3 1-2248 0-7492 1 -2948 0-9827 
0-4 1-3121 0-7803 1-4189 1-0472 
0-5 1-4045 0-8091 1-5549 1-1099 
0-6 1-5005 0-8342 1-7003 1-1672 
0-7 1-5983 0-8547 1-8508 1-2154 
0-8 1-6957 0-8696 2-0079 1-2510 
0-9 1-7907 0-8785 2-1450 1 -2722 
1-0 1-8814 0-8814 2-2789 1-2789 
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TABLE IV 
SUMMARY OF PLAN FORM CHARACTERISTICS 
dC 
S tan 9, tan A tan Ltan ¢, x K 
a 


Type A: Leading edge straight for x < 1:0 


0-5 0-876 0-677 2-094 2-605 0-830 1-009 
0-7 1-622 1-119 3-088 3-392 0-966 1-013 
1-0 4-040 2-265 5-079 4-353 1-360 1-047 
Type B: Leading edge straight for x <= 7? 

°5 2:034 1-151 2-603 3-238 0-828 1-006 
0-7 2-688 1-519 3-434 3-773 1-007 1-016 
1-0 4-040 2-265 5-079 4-353 1-360 1-047 

Type C: Fully curved leading edge 

0°5 1-332 0-881 2-334 2:904 | 0-793 1-012 
0-7 2-034 1-279 3-218 3:540 | 0-978 1-013 
1-0 4-049 2-265 5-079 4-353 1-360 1-047 


been taken as 45° for simplicity, but the corresponding shapes for other values of 
can be obtained at once by multiplying all lateral (vy) dimensions by cot ¢, while 
leaving all longitudinal (x) dimensions unchanged. For each value of m, it is clear 
that the aspect ratio A and the ratio of semi-span to root chord (s) are proportional 
to cot ¢,; the corresponding values of A tan¢, and stan, are given in Table IV 
(Type A) and are shown in Fig. 6 (Type A) plotted against 7 ,*=1/m,. 


2.4.2. Type B 

A more general class of leading-edge shapes can be obtained by starting (as in 
Section 2.4.1) with a straight portion OP,, x=m,», but now taking the point P,, out- 
board of which the singularity is to be held constant to be 


If x, <1 (as in all examples considered here), the shape of the curve P,B (Fig. 2) 
required to keep »=constant=o, for x,<x<1 may easily be found from 
equation (16) by setting G= H (1 —k*)=1 (compare the two forms of equation (8)) 
This gives 


where K=(x-—1)/» as before (k is, of course, negative). The value of «, may be 
found from equation (8), which gives 


No 21) 
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(the value of » when x= 1) can then be found by substituting k =k, =(x, —1)/n, in 
equation (20) as follows: 


n= V(l+e,7(k, exp | (22) 


Equations (20), (21) and (22) determine the shape of the curve PB (x, < x < 1), and 
the remainder of the leading edge (x > 1) can then be found as before. 


The leading-edge shapes just described depend on two basic parameters, m, (the 
taper of the inner part of the plan form) and x, (the position at which the leading edge 
becomes curved). To provide a contrast to the plan forms of Type A described in 
Section 2.4.1, for which x, was taken as 1 and m, varied, we have chosen as an 
illustration the second family of curves for which m, = 1 (i.e. there is no initial taper) 
and x, is varied. In this case x,=7,=«,°. Three examples, having x,=0:5, 0-7 and 
1-0, are given in Table III and shown in Fig. 5 (Type B); the last example is identical 
with the corresponding case in Section 2.4.1 (Type A). The corresponding values of 
A tan ¢, and s tan 4, are given in Table IV and plotted against «,* (=x,) in Fig. 6. 
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2.4.3. Type C 


The plan forms described in Sections 2.4.1 (Type A) and 2.4.2 (Type B) may 
suffer from certain disadvantages when A tan ¢, is small. For plan forms of Type A 
this implies that the initial taper parameter m, is large, and it is then found (see 
Section 3.2) that the local section lift coefficients increase rather rapidly over the 
inner part of the span. Plan forms of Type B are satisfactory in this respect, but 
may prove undesirable from the structural point of view since they have no initial 
taper, while the rapid increase of leading-edge curvature at P, may have aerodynamic 
disadvantages. These difficulties could be overcome by a proper choice of m, and 
x, (cf. Section 2.4.2) but, as an alternative, it is possible to prescribe the variation of 
o over the whole of the segment OB in such a way that it increases less rapidly in the 
spanwise direction than for the simple tapered plan forms of Type A (Section 2.4.1); 
this leads to leading-edge shapes which are curved over the entire span. 


The simplest expression for « which has the required property is 


Equation (8) then becomes 
x 
which has the solution 
(25) 
and . ‘ (27) 


In this way a third one-parameter family of leading-edge shapes (Type C) has been 
developed. Three examples, corresponding again to 7 ,?=0-5, 0-7 and 1-0, are given 
in Table III and shown in Fig. 5(C); the plan form with o,?=A=1 is again identical 
with the corresponding plan form of Type A (Section 2.4.1) or Type B (Section 2.4.2). 
The values of A tan¢, and stan, are given in Table IV and plotted against o,* 
in Fig. 6. 


2.5. HIGHER APPROXIMATIONS TO THE FUNCTION H (x) 


It will be recalled that the leading-edge shapes developed in Sections 2.3 and 
2.4 were based on an approximation (17) or (19) for the function H. When a plan 
form has been obtained in this way it is important to verify that the chosen 
approximation is in fact a good one and, if necessary, to modify the shape 
accordingly. This can be done by substituting the first approximation to H into the 
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right-hand side of equation (6); the left-hand side will then give the required second 
approximation to H. Only one case has in fact been investigated in this way, namely 
the plan form having m,=o,7=1. The results are shown in Fig. 7, in which the first 
approximation (19) to H and the second approximation obtained by numerical 
integration of equation (6) are plotted against the parameter k. The agreement 
between the two approximations is in general very good. The first approximation 
appears to overestimate the value of H; as expected, the difference between the two 
increases as k approaches unity (towards the wing tips), but it is less than 5 per cent 
when k=0:7 (y/s=0-90) and even when k=0-9 (y/s=0-99) it is only 13 per cent 
(of the second approximation). The actual errors in H should be less than these 
figures, since successive approximations will probably oscillate about the true value. 


In view of the fact that linearised potential theory is inevitably liable to error 
(due largely to viscous effects) in the immediate neighbourhood of the tips of 
sweptback wings, it is felt that the accuracy obtained with the present approxima- 
tion should be adequate for the purpose required. In addition, for plan forms of 
smaller aspect ratio, for which o, <1, the approximation should be still better, 
because the mean taper is greater and the terms neglected on the right-hand side of 
equation (17) will therefore, in general, be smaller. 


3. Wing Loading Distributions at M=1 


This section reviews the most important theoretical characteristics of the wings 
designed according to the methods of Sections 2.3 and 2.4. 
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3.1. CHORDWISE LOADING CURVES 


For a plane wing of given plan form, the wing loading coefficient / (the 
difference between the pressure coefficients on the upper and lower surfaces) can be 
found at any point by means of equations (3) (for x < 1) and (4) (for x > 1). For 
plan forms conforming to the design method of Section 2.3, the “unknown” function 
H in equation (4) has been shown to be given with sufficient accuracy by equation 
(19); H depends only on k, the ratio of the ordinate of the trailing edge to that of the 
leading edge at the particular value of x under consideration, and on the parameter 
«,, which represents the strength of the leading-edge singularity and is a constant for 
any one of these plan forms. 


The loading distribution at several spanwise stations has been calculated in this 
way for all the examples mentioned in the previous section. Some results are shown 
in Fig. 8, which gives the chordwise loading curves for the three wings (Types A, B 
and C) having «,* = 4 and for the wing with «7, =1. The most noteworthy feature of 
all these results is that in the tip region (outboard of the station »=»,), the chordwise 
load distribution, when plotted against £=(x—.x,)/c, is almost independent of span- 
wise position. This phenomenon could in fact have been anticipated; it is known 
(see Ref. 1) that, for plane wings in supersonic flow, the load distribution outboard 
of the station at which the Mach line from the root trailing edge crosses the leading 
edge is similar to the “flat-plate” loading ] oc /[(1 —&)/€] (cf. Fig. 8) and can be 
obtained approximately from it provided that the factor of proportionality is known. 
This factor clearly depends only on the strength o of the inverse square root 
singularity at the leading edge so that, if « is held constant, then the loading curves 
will not alter appreciably across the span. Wings designed in this way should thus 
possess the desirable feature that the loading remains constant as the tip is 
approached along lines of constant chordwise position; this is illustrated in Fig. 5 


which shows the lines of constant loading for all the plan forms downstream of the 
line x=1. 


3.2. SPANWISE VARIATION OF LOCAL LiFT COEFFICIENT 


A similar feature is present in the curves of local lift coefficient against spanwise 
position, shown in Fig. 9. In all the examples considered, the lift coefficient is nearly 
constant over the outer half of the wing span; in most cases there is in fact a very 
small decrease in C, as the tip is approached. There is, however, a noticeable 
difference in the behaviour of the three types of plan form over the inner part of the 
span. The wings which have straight tapered leading edges for x <1 (Type A. 
Section 2.4.1) show an appreciable increase in C,, in the spanwise direction, and for 
the smaller values of «,* (or aspect ratio), when the initial taper is large, this increase 
is sufficiently marked to form an undesirable feature unless other steps are taken to 
counteract it. The wings which are initially untapered (Type B, Section 2.4.2), on 
the other hand, all show a definite decrease in C,, (after a very small initial increase), 
the overall variation in C,, is less than for a corresponding wing of Type A. The 
third type of plan form (Type C, Section 2.4.3), which has a completely curved 
leading edge and is, in all respects, intermediate between Types A and B, turns out 
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to be quite a good approximation to a “constant C,”” plan form (cf. Brebner ”’); the 
overall variation in C,, is in no case greater than +5 per cent of the mean value. 


3.3. SPANWISE LOADING AND INDUCED DRAG 


The corresponding spanwise loading curves are shown in Fig. 10, in which the 
ratio 1’=(2sc/S) (C,./C,) (where C,, is the overall lift coefficient) is plotted against 
spanwise position y/s. It can be seen that the shape of the loading curve depends 
very little on the type of plan form chosen. For the smaller aspect ratios (with 
o,* = 4) the loading is nearly elliptic, and even for the highest aspect ratio considered 
(A tan ¢, = 5-08), where the loading at the root is about 12 per cent less than its 
maximum value, the effect on the vortex drag is not serious. This can be seen from 
Fig. 11, where the values of the vortex drag factor K, defined by 


Cy 
2. =A (28) 
C, 
are plotted against A tan ¢, for all the examples considered. These factors have been 
calculated from the loading curves by the standard method of Fournier analysis. The 
results show a certain amount of scatter, presumably due to errors in the numerical 
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evaluation of the local lift coefficients, but the resulting errors in K should not 
exceed } per cent. Within this accuracy there is no significant difference between the 
three types of plan form, and all the results may be approximated by the common 
formula 


K=1+0-0018(A tan ¢,)’. j (29) 


Thus, over a wide range of aspect ratios of practical interest (0< A tan ¢, <5), 
the vortex drag factor does not exceed 1-05. It seems therefore that all the plan 
forms developed here are good approximations to theoretical ‘“‘minimum vortex 
drag” plan forms, particularly for the smaller values of A tan¢,, and it would seem 
pointless to seek any closer approximations to this ideal. A much more important 
question is, whether these low values of K would be obtained in practice. This is a 
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point which can only be settled by experiment, but it seems probable that wings 
which are designed, either by choice of plan form or other methods, to have constant 
pressure along lines of constant chordwise position will be more satisfactory in this 
respect than more conventional sweptback wings, where the high suction peaks near 
the wing tips are a major factor in causing early boundary layer thickening and 
separation in the tip regions. y 


3.4. OVERALL LIFT AND CENTRE OF PRESSURE 


The overall lift coefficient C,, may be obtained either by numerical integration of 
the spanwise loading curves or, more accurately, in the following way. 


Mangler” has shown that the integral of the local loading coefficient / over a 
strip in the spanwise direction is 


Jldy 4z2Hy, (1 —«) (30) 
dx 
(note that H (1 —«)=1 for x <= 1). 
IC, 
Hence k)ydy , 
da 
ome 1 
dC, ( dn) 
an ¢ dk | 32) 


This can be written in the form 


da A ( [Ha k) Guo? pak | (33) 


0 


where n,/», is given by equation (16) or Table II. Since both »,/», and s/y, depend 
only on o,, it follows that (dC, /dz) (471A) is also a function of «, only. Further, for 
any wing of a given type, «, determines the value of A tan ,, so that (dC,/ dz) tan 9, 
is a function of A tan ¢,, as indeed could be deduced from similarity considerations. 


The values of (dC,,/dz) tan, are given in Table IV and are plotted against 

A tan, in Fig. 12, while the common values of (dC,,/dz)/(47A) are plotted against 
o, in Fig. 13. In both these figures the left hand (small aspect ratio) portions of the 
curves have been estimated in the following way. In the limit as 7,’ (or A tan 9,) 
approaches zero, it may be shown that the integral in equation (33) approaches the 
value 

1 

Ja-K) ak ~ 030. 
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The corresponding spanwise loading curves are shown in Fig. 10, in which the 
ratio 1’=(2se/S) (C,./C,) (where C, is the overall lift coefficient) is plotted against 
spanwise position y/s. It can be seen that the shape of the loading curve depends 
very little on the type of plan form chosen. For the smaller aspect ratios (with 
o ,?=4) the loading is nearly elliptic, and even for the highest aspect ratio considered 
(A tan ¢, = 5-08), where the loading at the root is about 12 per cent less than its 
maximum value, the effect on the vortex drag is not serious. This can be seen from 
Fig. 11, where the values of the vortex drag factor K, defined by 
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are plotted against A tan ¢, for all the examples considered. These factors have been 
calculated from the loading curves by the standard method of Fournier analysis. The 
results show a certain amount of scatter, presumably due to errors in the numerical 
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evaluation of the local lift coefficients, but the resulting errors in K should not 
exceed } per cent. Within this accuracy there is no significant difference between the 
three types of plan form, and all the results may be approximated by the common 
formula 

K=1+0-0018 (A tan ¢,)’. (29) 


Thus, over a wide range of aspect ratios of practical interest (0< A tan, <5), 
the vortex drag factor does not exceed 1-05. It seems therefore that all the plan 
forms developed here are good approximations to theoretical ‘‘minimum vortex 
drag” plan forms, particularly for the smaller values of A tan, and it would seem 
pointless to seek any closer approximations to this ideal. A much more important 
question is, whether these low values of K would be obtained in practice. This is a 
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point which can only be settled by experiment, but it seems probable that wings 
which are designed, either by choice of plan form or other methods, to have constant 
pressure along lines of constant chordwise position will be more satisfactory in this 
respect than more conventional sweptback wings, where the high suction peaks near 
the wing tips are a major factor in causing early boundary layer thickening and 
separation in the tip regions. 


3.4. OVERALL LIFT AND CENTRE OF PRESSURE 


The overall lift coefficient C;, may be obtained either by numerical integration of 
the spanwise loading curves or, more accurately, in the following way. 


Mangler®’ has shown that the integral of the local loading coefficient / over a 
strip in the spanwise direction is 


dx 


(note that H (1 —«)=1 for x <= 1). 


s 


dC, Ly, 


Hence 4 H(1—«)ydy . (31) 


or tan @,- tk | (32) 


This can be written in the form 


1 
where »,/», is given by equation (16) or Table II. Since both »,/», and s/y, depend 
only on o,, it follows that (dC,,/dz)/(474A) is also a function of «, only. Further, for 
any wing of a given type, «, determines the value of A tan 4, so that (dC,,/ dz) tan 9, 
is a function of A tan ¢,, as indeed could be deduced from similarity considerations. 


The values of (dC,,/dz) tan, are given in Table IV and are plotted against 
A tan 4, in Fig. 12, while the common values of (dC,,/dz)/(47A) are plotted against 
o, in Fig. 13. In both these figures the left hand (small aspect ratio) portions of the 
curves have been estimated in the following way. In the limit as 7,* (or A tan 9,) 
approaches zero, it may be shown that the integral in equation (33) approaches the 
value 
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Also, as o,” approaches zero, ) ~ 1-o,?, 
\ 
so that equation (34) takes the approximate form 


dC, 


The relationship between «,’* and the aspect ratio A, when both are smal!, may 
be found by consideration of the limiting shapes taken by the three types of plan 
form. These are, respectively, a standard delta (Type A), for which 

A tang, ~ 


a half ellipse (Type B), for which 
8 
A tan¢, > — . (36) 


and a parabolic “Gothic” delta (Type C), for which 


3 
A tan ¢, = . (37) 


V 


Substituting these results into equation (34), we obtain the following approximations, 
when A tan 4, is small, for the lift curve slopes of the three types of wings : 


dC, 


Tvpe A: ig 4xA (1—0°10 A tan¢,) (38) 

Type B: ~ (1 —0-062 4? tan? (39) 

Type C: (1—0°067 A’ tan’ . . (40) 


These expressions appear to be valid for A tan¢, < I. 


The position of the centre of pressure can be found in a similar way. Measuring 
the pitching moment about the apex (x =0) and taking the root chord as the standard 
length, we have 


1 


da (stang,)?L) xa — dk dk |. 


Then x= —(dCy/dC;) gives the theoretical distance of the centre of pressure (or 
aerodynamic centre) downstream from the apex as a fraction of the root chord. 
Values of X are given in Table IV and are plotted against A tan ¢, in Fig. 14. The 
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values shown for A tan #,=0 depend only on the limiting shapes to which the three 
types of plan form tend as the aspect ratio or trailing-edge sweep approaches zero. 
Thus for wings of Type A (standard delta) the limiting value of x is 2/3; for Type B 
(half ellipse) it is 1/3; and for Type C (parabolic Gothic delta) it is 7/15. 


4. Discussion 


The method suggested in the present paper enables a wing plan form to be 
designed to have any desired variation in the strength of the singularity in the load 
distribution along the leading edge, for a flat wing at M=1. Outboard of the point 
at which the Mach line from the root trailing edge meets the leading edge, the load 
distribution in a chordwise direction is everywhere similar to a two-dimensional “flat 
plate” distribution, so that spanwise control of the leading-edge singularity also 
controls the loading behind it. In particular, when the leading-edge singularity is 
held constant the chordwise loading also becomes independent of spanwise position. 
This property has been used in the design of three families of wings, in each of 
which the aspect ratio is varied by altering the value of the constant singularity 
Strength in the tip region. The three families differ chiefly in their shape and 
aerodynamic behaviour over the inner part of the span. 


Wings of Type A (Section 2.4.1) have straight tapered leading edges up to the 
point (x= 1) at which the Mach line from the root trailing edge meets the leading 
edge; they show a tendency for the local lift coefficient to increase over the inner part 
of the span, which may be undesirable if it becomes too large, as is the case when the 
initial taper is high and the aspect ratio (or more precisely A tan 9,) is low. Wings 
of Type B (Section 2.4.2) are initially untapered, and the leading-edge singularity is 
held constant from a point P, that is farther forward than for Type A. As a result, 
the load distribution is more favourable and in particular the local lift coefficient 
decreases slightly over most of the span but, for the smaller aspect ratios, the plan 
forms show a marked change in leading-edge curvature at P,. and this may be 
undesirable at low speeds and high incidences because of the tendency for a vortex 
to form at P,. Wings of Type C (Section 2.4.3), which have fully curved leading 
edges, are in every sense intermediate between those of Types A and B and may be 
preferable in some respects, apart from the additional geometrical complication of 
the plan form shape. 


The differences noted between the three types of plan form are only significant 
when the aspect ratios are small, and in fact they disappear altogether in the case of 
the highest aspect ratio here considered (A tan ¢,=5-08). In Section 2.4.2 a more 
general type of plan form, combining the properties of Types A and B, was also 
discussed, though not investigated in detail. Here the leading edge is initially 
straight and tapered (as for Type A), while the singularity is held constant from an 
arbitrary point P, (as for Type B); such a plan form might be useful in a case where 
one of the simpler types was found unsuitable. 


All the calculations described here have been for plane wings alone. The 
question naturally arises: What will be the effect of mounting these wings on a 
fuselage? When the fuselage is slender and of circular cross section, Mirels® has 
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shown how to find the interference effects by means of a Joukowski transformation 
in the cross-flow planes; this reduces the problem to that of calculating the loading 
on a wing alone of distorted plan form. Although this is a comparatively straight- 
forward matter, it is not possible to use any of the results of the present paper 
directly; but it is clear that, if the fuselage radius is small compared with the wing 

semi-span, then the load distribution in the tip regions will not be greatly altered by 
the presence of the fuselage, so that the principal features of the plan forms designed 
here should not be seriously affected. 


The entire design method is limited by the use of linearised lifting-surface 
theory, which neglects any interaction between thickness and incidence and ignores 
the effect of viscosity. Apart from these limitations there is no particular reason to 
suspect the use of the simplified (sonic) form of the theory, particularly if the wing 
design is successful in eliminating some of the shock waves which might otherwise 
occur. The whole matter is clearly in need of experimental verification, and a 
programme of tests with this object, among others, in view is in progress at the N.P.L. 
The wing chosen has the plan form described in Section 2.4.1, and shown in Fig. 5, 
having o,?7=1, A tan¢,=5-08, with ¢,=55°, and will be tested with and without a 
fuselage over a Mach number range from 0°8 to 1:2. As well as providing a useful 
experimental check on the correctness of the theory and design method, these tests 
should also help to determine how far on either side of the design Mach number 
(M =1-0) it is possible to extend the basic theoretical results. 
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Turbulent Flows and Heat Transfer. Volume V of High Speed Aerodynamics and Jet 
Propulsion. Edited by C. C. Lin. Oxford University Press, London, 1960. 550 pp. 
Illustrated. £5 5s. Od. 


It is difficult to avoid the irreverent reflection, when reading this volume, that the 
Princeton series on High Speed Aerodynamics and Jet Propulsion has got slightly out 
of hand. One of the title pages records the toll of general editors during the 12 years 
since the project was begun, and some of the changes of plan are apparent from a com- 
parison of the two lists of titles and editors of the separate volumes, one on a title page 
and one on the dust-jacket. The majestic comprehensiveness of the series, which must 
have seemed a noble and inspiring goal in 1949, now begins to look a little like a mill- 
stone, as contributions written promptly become obscelescent while waiting on slower 
contributors to the same volume and as the master-plan suffers more modification. 
Four more volumes remain to be published. 


The present volume has not escaped the penalties associated with the series. It is 
a heterogeneous mixture of contributions on a variety of topics, with very little unity 
of purpose, subject matter or style; some of the articles look a little out of touch with 
current work, and there are some omissions which, for a book published in 1959, 
seem to be serious. (For instance, the important and interesting problem of generation 
of noise by turbulent flows—the inclusion of which would have helped to give the book 
contact with “high speed aerodynamics and jet propulsion ”’—is hardly mentioned.) 
Some of the contributions are useful surveys (especially useful are those that gather 
together and collate large numbers of experimental observations), and these undoubtedly 
give the book value; however, they vary considerably in scope and quality and a thorough 
review would need to analyse them one by one in more detail than is appropriate here. 
All that I propose to do is to list the titles and authors of the separate articles, and to 
point out to a prospective purchaser that these articles have the general character to 
be expected of a brief review of the named subject by the named author, that they can 
be read as entirely self-contained reviews, and that the whole volume does not provide 
a complete picture of anything in particular. 


The nine sections of the book are as follows : — 


Transition from Laminar to Turbulent Flow, by H. L. Dryden (72 pages) 
Turbulent Flow, by G. L. Schubauer and C. M. Tchen (121 pages) 
Statistical Theories of Turbulence, by C. C. Lin (58 pages) 

Conduction of Heat, by M. Yachter and E. Mayer (34 pages) 


Convective Heat Transfer and Friction in Flow of Liquids, by R. G. Deissler and 
R. H. Sabersky (51 pages) 


F. Convective Heat Transfer in Gases, by E. R. van Driest (89 pages) 
G. Cooling by Protective Fluid Films, by S. W. Yuan (61 pages) 
H. Physical Basis of Thermal Radiation, by S. S. Penner (13 pages) 
I. Engineering Calculations of Radiant Heat Exchange, by H. C. Hottel (39 pages) 
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AERO POWER 


LIFT and THRUST from ONE ENGINE-- 
the revolutionary BS 53 turbofan 


The new Bristol Siddeley BS 53 high-ratio turbofan 
represents a major engineering breakthrough in the 
field of aircraft propulsion. Because the BS 53 can be 
fitted with movable nozzles evenly disposed around 
the centre of gravity which give directional control 
to its total thrust, this remarkably versatile engine 
provides the airframe designer with a single power 
source for all conditions of flight. 


VTOL, STOL and normal take-off capability in the same 
aircraft. Suitable for subsonic or supersonic applications. 


High cold-flow ratio gives high thrust for low weight, 
low fuel consumption and low noise level. 


Conventional installation. Operational simplicity. Indee 
pendent of all fixed ground installations. 


Problems of ground erosion, recirculation and ingestion 
of debris reduced to a minimum. 

The Bristol Siddeley BS 53 has already been selected 
for, and is now undergoing trials in the Hawker P 1127 
VTOL light strike fighter and is supported by the 
Mutual Weapons Development Programme for NATO. 


BRISTOL SIDDELEY ENGINES LIMITED 


TURBOJETS—TURBOFANS—TURBOPROPS——-RAMJETS—ROCKET ENGINES—MARINE AND INDUSTRIAL GAS TURBINES— 
MARINE, RAIL AND INDUSTRIAL DIESEL ENGINES—PISTON ENGINES—PRECISION ENGINEERING PRODUCTS. 
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NOTICE TO CONTRIBUTORS 


1. All communications should be addressed to the Editor, The Aeronautical 
Quarterly, 4 Hamilton Place, London, W.1. 


2. Papers for consideration should be typewritten with double spacing and be 
accompanied by a summary of not more than 250-300 words. Papers must be 
as short as possible, otherwise considerable delay in publication will result, and 
in any case they should not be longer than twenty foolscap pages. Full 
descriptions and complete results of work, as written for an internal report, are 
not suitable for publication in The Aeronautical Quarterly, and such reports 
need to be rewritten in concise form. Only necessary equations should be 
included in the text; routine mathematics should be omitted or given as an 
Appendix to the paper. [Illustrations should be reduced to a minimum. The 
titles of papers should be kept short—with sub-titles if necessary. 


3. Only simple symbols and formulae should be typewritten. All others should 
be written clearly in ink. Greek letters should be indicated and a clear distinction 
drawn between capital and small letters. Particular attention should be paid 
to the use of alpha and a, kappa and k, mu and wu, nu and v, eta and nm. The 
letter / must be distinguished from the numeral 1, and the letter O from the 
zero 0. 

The use of dots, bars, and so on, over letters, or the use of dots as multi- 
plication signs and bars for brackets should be avoided as far as possible. 

Suffixes and prefixes must be indicated clearly and complicated suffixes 
and prefixes should be avoided. 


4. A list of all symbols used, both in the illustrations and text, must be given 
whether they are standard or not. 


5. Clear diagrams or sketches must be submitted for redrawing by the Society’s 
draughtsman. 
Photographs must be not less than half plate size and must be clear black 
and white glossy prints. 


6. References should be given in the form: — 
Warp, G. N. Calculation of Downwash Behind a Supersonic Wing. The 
Aeronautical Quarterly, Vol. I, p. 35, May 1949. 


7. The copyright of every paper printed in The Aeronautical Quarterly is the 
property of The Royal Aeronautical Society. Permission to reprint or to use 
any paper will not be unreasonably refused. 


8. It is a condition of publication that the author shall have obtained consent where 
necessary to use any material in his paper which is copyright or the property 
of any other person or his employers. 


9. Authors will be entitled to 20 reprints of their papers free of charge. This 
number will be shared between joint authors. 


The Aeronautical Qua is published at the Offices of The Royal Aeronautical Society, 
4 Hamilton Place, Piccadilly, W.1, at 15s. Od. net for a single number to non-members 
and 7s. 6d. net for a single number to members of The Royal Aeronautical Society. Instructions 
for Contributors are given above. None of the papers must be taken as expressing the 
opinion of the Council of The Royal Aeronautical Society unless such is definitely stated to 
be the case. All communications for publication or advertisement in The Quarterly should be 
addressed to—The Editor, The Aeronautical Quarterly, 4 Hamilton Place, W.1. 


Telephone: Grosvenor 3515. Telegraphic Address: Didaskalos, Audley, London. 
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THE AERONAUTICAL QUARTERLY 


It is the aim of The Aeronautical Quarterly to attract not only original papers contributing 
to aeronautical science and engineering, and papers developing new or improved methods of 
analysis and experimental techniques, but also papers on allied sciences which have a bearing 
on aeronautical problems. The Aeronautical Quarterly is open to authors of any nationality 
and is not restricted to members of the Society. 


Most papers can be published within nine months of receipt. Short papers, of not more 


than 2,000 words, are published much more quickly. 


A selection of the more specialised papers read before the Society during the current 


lecture session will be published in The Aeronautical Quarterly. 


VOLUME XII 


Volume XII of The Aeronautical Quarterly will be published in four parts, in February, 
May, August and November 1961. The prices, including postage and packing, are: — 


To Members of the Society To Non-Members 


£s. 
Per Part 8 0 ($1.15) 15 6 ($2.20) 


Subscription (4 parts) £220 ($4.50) 3 2 0 ($8.70) 


PRINTED BY THE LEWES PRESS WIGHTMAN & CO. LTD., LEWES, SUSSEX, ENGLAND, AND PUBLISHED 
BY THE ROYAL AERONAUTIC..L SOCIETY, 4 HAMILTON PLACE, LONDON W.1, ENGLAND. 
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